{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Fit CoFactor model to the binarized ML20M"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "import itertools\n",
    "import glob\n",
    "import os\n",
    "import sys\n",
    "os.environ['OPENBLAS_NUM_THREADS'] = '1'\n",
    "\n",
    "import numpy as np\n",
    "import matplotlib\n",
    "matplotlib.use('Agg')\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "\n",
    "import pandas as pd\n",
    "from scipy import sparse\n",
    "import seaborn as sns\n",
    "sns.set(context=\"paper\", font_scale=1.5, rc={\"lines.linewidth\": 2}, font='DejaVu Serif')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import cofacto\n",
    "import rec_eval"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Construct the positive pairwise mutual information (PPMI) matrix"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Change this to wherever you saved the pre-processed data following [this notebook](./preprocess_ML20M.ipynb)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "DATA_DIR = '/hdd2/dawen/data/ml-20m/pro/'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "unique_uid = list()\n",
    "with open(os.path.join(DATA_DIR, 'unique_uid.txt'), 'r') as f:\n",
    "    for line in f:\n",
    "        unique_uid.append(line.strip())\n",
    "    \n",
    "unique_sid = list()\n",
    "with open(os.path.join(DATA_DIR, 'unique_sid.txt'), 'r') as f:\n",
    "    for line in f:\n",
    "        unique_sid.append(line.strip())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "111148 11711\n"
     ]
    }
   ],
   "source": [
    "n_items = len(unique_sid)\n",
    "n_users = len(unique_uid)\n",
    "\n",
    "print n_users, n_items"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def load_data(csv_file, shape=(n_users, n_items)):\n",
    "    tp = pd.read_csv(csv_file)\n",
    "    timestamps, rows, cols = np.array(tp['timestamp']), np.array(tp['uid']), np.array(tp['sid'])\n",
    "    seq = np.concatenate((rows[:, None], cols[:, None], np.ones((rows.size, 1), dtype='int'), timestamps[:, None]), axis=1)\n",
    "    data = sparse.csr_matrix((np.ones_like(rows), (rows, cols)), dtype=np.int16, shape=shape)\n",
    "    return data, seq"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "train_data, train_raw = load_data(os.path.join(DATA_DIR, 'train.csv'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "watches_per_movie = np.asarray(train_data.astype('int64').sum(axis=0)).ravel()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The mean (median) watches per movie is 597 (48)\n"
     ]
    }
   ],
   "source": [
    "print(\"The mean (median) watches per movie is %d (%d)\" % (watches_per_movie.mean(), np.median(watches_per_movie)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "user_activity = np.asarray(train_data.sum(axis=1)).ravel()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The mean (median) movies each user wathced is 62 (33)\n"
     ]
    }
   ],
   "source": [
    "print(\"The mean (median) movies each user wathced is %d (%d)\" % (user_activity.mean(), np.median(user_activity)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "vad_data, vad_raw = load_data(os.path.join(DATA_DIR, 'validation.csv'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAbQAAAE5CAYAAAD4GtVAAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3XdYVNf2N/Dv0KUjiAq2WCKCCImgUWNQrPFiJDEq1thN\nLCGIBTUieGNvKMTYMaKCBQHFxBIQ9XoTkVwbP0QFOyBlkCY6QznvH7ycOMLAYZgZZuasz/P4PM5p\nsxaDrtn77LO3gGEYBoQQQoia02rqAAghhBB5oIJGCCFEI1BBI4QQohGooBFCCNEIVNAIIYRoBCpo\nhBBCNAIVNEIIIRqhUQXt+fPn8oqDEEIIaZQGFbS8vDxkZmayfxYvXqyouAghhJAG0eFyUGxsLAID\nA1FSUsJuYxgGAoFAYYERQgghDcGpoG3fvh0rV65E9+7doa+vD6CqoPn6+io0OEIIIYQrTgXNxsYG\nX3zxRY3tQUFBcg+IEEIIkQWne2i9evXCw4cPa2zftWuX3AMihBBCZCHgMtv+5MmTkZycjPbt28PU\n1JTdnpqaisTERIUGSAghhHDBqcsxLy8P/v7+EtsYhsHevXsVEhQhhBDSUJwK2uzZs/Hll1/W2G5k\nZCT3gAghhBBZcOpyrHb37l1kZmbCxsYGjo6ODXqj8vJyHD58GHFxcQCAsrIyeHt7o0+fPgCATz/9\nFB07dmQfB+jZsye8vb3Z8w8cOIDY2FgIBAKMGDECM2bMYPcVFhbCz88PhYWFqKyshL+/P+zt7RsU\nHyGEEDXHcJCXl8eMHTuW6dq1K/tnzJgxTG5uLpfTGYZhmBcvXjCDBg1iSkpKGIZhmGvXrjHOzs5M\ndnY2wzAM4+fnJ/Xcy5cvM59//jkjFosZkUjEDB8+nElISGD3+/j4MMHBwQzDMMxff/3FuLm5MWKx\nmHNshBBC1B+nUY5r165Fz5498ccff+DWrVu4ePEievbsibVr13IunEZGRvj+++/Zbsq+fftCX18f\nN2/erPfc48ePw8PDA7q6utDT08PIkSMREREBACgoKMC5c+fw9ddfAwB69+4NXV1dJCQkcI6NEEKI\n+uNU0LKzs7FkyRK0adMGBgYGaNu2LZYuXYrs7GzOb2Rubl7jWbaysjJYWloCANLT0zF79mxMmDAB\n/v7+ePXqFXvc3bt30bFjR/Z1586dkZycDABISUmBvr4+WrVqxe7v1KkTu58QQgg/cCpo5eXltW6v\nqKiQ+Y2vX78OW1tbuLi4AAC6dOmCoKAgHD16FMbGxpg1axZ7rFAohImJCfvaxMQE+fn57D5jY2OJ\na5uamkIoFMocGyGEEPXDqaDZ2dlh4cKFuHHjBtLT05GYmIiFCxeia9euMr2pSCTCtm3bsGHDBnbb\nmjVrYGhoCABYsGABUlJScOfOHXZ/XfNG1raP4T7WhRBCiAbgNGx/8eLFWLZsGSZPngyBQACGYTB0\n6FAsWbJEpjf19/fH9OnT0a1bt1r3N2vWDGZmZsjMzESPHj3QvHlzFBUVsfuLi4thYWEBALC0tERx\ncbHE+UVFRfjwww+lvv+Pv1zD0N7t0cextUzxqxttbQEqKvhR4PmUK8CvfPmUK8CvfHV1teVyHU4F\nzcjICDt27IBQKERGRgZsbW3Ze18NtX79ejg5OWHo0KEQi8UQCoV4+vQpjI2N0b17dwBV99aKiorQ\nsmVLAICjoyMeP37MXiMtLY19bMDe3h5isRjZ2dns8enp6Rg9erTUGG6n5eF2Wh7menaHi521THmo\nE3NzQxQUlDZ1GErBp1wBfuXLp1wBfuXbooVJ/Qdx0KD10CwtLdGjRw+Zi9mePXtQXl4OT09PlJaW\n4tmzZ4iMjERWVhbCw8PZ48LCwtCuXTv06NEDAODl5YWzZ89CLBZDJBIhNjYW48ePB1A12GT48OE4\nefIkACAxMRHl5eVwc3OrN56EWxky5UEIIUT1cGqhycOTJ0+wdetWCAQCHDlyhH2Aet68eejTpw9u\n3LiBiRMnorKyEsbGxti1axe0tauaof3790d6ejq8vLwgEAgwZswYfPbZZ+y1/f394efnhwkTJqCy\nshK//PIL9PT0lJUaIYQQFdCgmUI0xUjfGACgLkcNxKdcAX7ly6dcAX7lK68uR6W10FRJj85W+LR7\nK14UM0II4YsG3UN7X2RkpLziUKqAmZ9QMSOEEA0jtYUWHR1d78n79++vczQhIYQQoixSC9qaNWtg\nZ2fHvr516xbat28PCwsLvHr1Co8ePcLHH3+slCAJIYSQ+kgtaH379sX27dsBANu3b4ePj49EAUtK\nSsKFCxcUHyEhhBDCgdR7aNXFDKhqnb3fGnNxccH9+/cVFxkhhBDSAJwGheTm5taY7Dc3Nxd5eXkK\nCYoQQghpKE7D9keOHIl//etfGDhwICwtLSEUCnHp0iXMnDlT0fERQgghnHAqaHPmzEHLli1x+vRp\n3L59G61bt8aPP/4IDw8PRcdHCCGEcML5wWpPT094enoqMhZCCCFEZpwfrL5y5Qrmz5+POXPmoKio\nCDt27GjUAp+EEEKIPHEqaCdOnIC/vz9atGiBnJwcGBsbo7KyEhs3blR0fIQQQggnnApaVFQUTp8+\njVWrVsHY2BhaWlr44Ycf8H//93+Kjo8QQgjhhFNB09bWhqmpaY3tAoFA7gERQgghsuB8Dy0+Ph7A\nP0Xs+vXr0NJq1NzGhBBCiNxwGuW4ePFiTJ06Febm5iguLsa//vUv5OTkIDQ0VNHxEUIIIZxwKmg9\nevTAuXPnEBMTg8zMTNjY2GDUqFGwtqYlWAghhKgGTgUtOjoanp6emDVrlsR2X19fbNmyRSGBEUII\nIQ3B6SbY4cOHa2w7e/YsEhIS5B0PIYQQIhNOBe3Zs2cIDg4GAJSWlsLPzw8BAQEwMjJSaHCEEEII\nV5y6HOfOnQsdHR2sW7cOCQkJsLS0RExMDK5fv67o+AghhBBOOBW0qVOnAgCCg4PRvn177N69GwKB\nAKWlpYqMjRBCCOFMakELCQmpsU0gEODhw4fYtGkTDA0NERUVhYkTJyo0QEIIIYQLqQUtIiIC/fv3\nr7H9k08+watXr/Dq1SuIRCKFBkcIIYRwJbWgDR8+HD/++GOdJwcFBck9IEIIIUQWAoZhmIae9Pbt\nW+jo6EBHh/NyaiqlrKwCBQX8uf9nbm7Im3z5lCvAr3z5lCvAr3xbtDCRy3U4Ddv/+eef0bt3b9y4\ncQNA1TB+Dw8Pmm2fEEKIyuBU0C5fvoyzZ8/C1dUVAPDhhx9i7969tB4aIYQQlcGpoBkYGMDKykpi\nW9u2bSFDbyUhhBCiEJwKmkgkQk5OjsS27OxsGuVICCFEZXAa1TF+/HiMHDkS7u7usLS0hFAoRHx8\nPFasWKHo+AghhBBOOBU0T09PmJiYICIiAnfu3IGNjQ3Wr1+PgQMHKjo+QgghhBPO4+4HDRqEQYMG\nSWxLTU2FnZ2d3INSdUmpOUi4lQEAGOBsCxc7WheOEEKaGqeClpmZWev2gIAAREREyDUgVZeUmoOd\n0cns65QnrzDXszsVNUIIaWJSC5qHhweOHj0KU1NTuLu7QyAQ1BjVKBAIFB6gqqlumb2/jQoaIYQ0\nLakFbfXq1TA1NQUAuLq6IiwsrMYxkydPVlxkGoK6JwkhRDmkDtv/+OOP2b+vW7eu1mPWr18v/4hU\n3ABnW07bgH+6J1OevELKk1fYGZ2MpNScWo8lhBDSOJyeQ2vTpk2t25ctWybXYNSBi5015np2h30H\nC9h3sKjz/pm07klCCCHyJ7XL8f0Rje9jGAZCoZDzG5WXl+Pw4cOIi4sDAJSVlcHb2xt9+vQBANy7\ndw+BgYHQ0tKCmZkZ1q9fDzMzM/b8AwcOIDY2FgKBACNGjMCMGTPYfYWFhfDz80NhYSEqKyvh7+8P\ne3t7zrE1lIudNXUdEkKIipFa0ExMTLB8+XKpJzIMI7UrsjbZ2dk4fPgwYmJiYGRkhP/+97+YO3cu\nzp8/DwsLC8ybNw8bNmyAq6srgoOD4e/vj+3btwMArly5gpMnTyImJgYMw2DUqFHo3Lkz3NzcAACB\ngYFwcHDA/Pnzcf36dcydOxcXL16Erq4u5/gUYYCzLVKevKqxjRBCiPxJLWg+Pj7o1atXnSf7+Phw\nfiMjIyN8//33MDIyAgD07dsX+vr6uHnzJnR0dKCtrc1OfjxmzBi4u7vj1atXsLCwwPHjx+Hh4cEW\nqJEjRyIiIgJubm4oKCjAuXPnEB8fDwDo3bs3dHV1kZCQgCFDhnCOTxGquydpUAghhCie1IJW3foB\ngBcvXiApKQmurq6wtbXFixcvkJaWhgEDBnB+I3Nzc3zxxRcS28rKytC8eXNcu3YNHTt2ZLe3atUK\nBgYGSElJQb9+/XD37l14eHiw+zt37ozw8HAAQEpKCvT19dGqVSt2f6dOnZCcnNzkBQ2ou3uSRkAS\nQoj8cBoUsm3bNly7dg1aWlWHN2vWDMePH0doaKjMb3z9+nXY2trC1dUVQqEQxsbGEvtNTU3Ze3RC\noRAmJv8sAGdiYoL8/Hx2X13nqioaAUkIIfLFqaBlZmZi06ZNaN26NQDA0tISISEhOH/+vExvKhKJ\nsG3bNmzYsIHdVttD2u8+yF3XQ9z1nauKaAQkIYTIF6epr3R0ah6mpaUl86ALf39/TJs2Dd26dQMA\nNG/eHKmpqRLHFBUVsWuwNW/eHEVFRey+4uJiWFhYAKgqrsXFxTXO/fDDD6W+v7a2AObmhjLFLi86\nOtq1blNEXNraWk2er7LwKVeAX/nyKVeAf/nKA6eCZmhoiPDwcIwdOxba2tqoqKjAsWPH0KxZswa/\n4fr16+Hk5IRhw4ZBLBZDKBTC0dERv//+O3tMVlYWRCIRHBwcAACOjo54/Pgxuz8tLQ2Ojo4AAHt7\ne4jFYmRnZ6Nly5YAgPT0dIwePVpqDBUVDAoKShscuzx92r0V7qTl1dimiLjMzQ2bPF9l4VOuAL/y\n5VOuAL/ybdHCpP6DOODU5bh8+XLs2bMHTk5O6NevH5ycnLBv374Gr4e2Z88elJeXw9PTE6WlpXj2\n7BkiIyPx2Wefoby8HElJSQCAyMhIDBkyBObm5gAALy8vnD17FmKxGCKRCLGxsRg/fjyAqsEmw4cP\nx8mTJwEAiYmJKC8vlxjUoooa8oA2IYSQ+gkYjjeb3rx5g0uXLiEzMxM2NjYYOHBgg1poT548wfDh\nw9n7XQzDQCAQYP78+Zg3bx77YLW2tjZMTU1rPFh98OBBnD59GgKBAB4eHpg2bRq77/0HqwMCAupc\n1qasrII333wAfn3T41OuAL/y5VOuAL/ylVcLjXNB0yRU0DQXn3IF+JUvn3IF+JWvUrscCSGEEFVH\nBY0QQohGoIJGCCFEI3AqaOvWrcPhw4cVHQshhBAiM04FLTY2Fr1791Z0LIQQQojMOBW0Hj16oEuX\nLjW2X758We4BEUIIIbLgVNC+/PJLHDlypMYUU/v371dIUIQQQkhDcZr66vvvvwcA/PTTT+y26gej\nCSGEEFXAqaA5OTlh69atEtsYhoGvr69CgiKEEEIailNBCwwMhK2tbY3tmzZtkntAhBBCiCw4FTQ7\nOzvk5+cjNjYWxcXFmD59Om7fvo1PPvlE0fERQgghnHAaFHLr1i0MGTIEx44dw9mzZyEQCPDLL7/g\nxIkTio6PEEII4YRTQdu6dStCQ0Nx9uxZWFpawsDAAHv37kV0dLSi4yOEEEI44VTQGIZBjx49AIAd\n2ainpwdt7ZqrLhNCCCFNgVNBE4vFePnypcS23NxciEQihQRFCCGENBSnQSHjx4/HyJEjMXjwYLx4\n8QKBgYH4448/sGTJEkXHRwghhHDCqaB5enrCxMQEERERaNasGV68eIHVq1dj4MCBio6PEEII4YRT\nQQOAQYMGYdCgQYqMhRBCCJEZp3tod+7cQUhICMrLy5GdnY2pU6di8uTJePTokaLjI4QQQjjhVNB2\n7NgBExMTaGlpISgoCAYGBnB1dZWY25HIT1JqDjZH3MTmiJtISs1p6nAIIUQtcOpyLCsrwzfffIPy\n8nJcunQJ586dg7m5OSZOnKjo+HgnKTUHO6OT2dcpT15hrmd3uNhZN2FUhBCi+ji10MrKygAAV69e\nhaOjI8zNzQEAurq6iouMpxJuZXDaRgghRBKnFtoHH3yAmTNn4v79+1i1ahXevn2Lo0ePQl9fX9Hx\nEUIIIZxwaqGtWrUKQ4YMwYoVKzB48GC8efMGpaWlmDdvnqLj450BzjVXNahtGyGEEEkChmEYWU9O\nTU2FnZ2dPONRirKyChQUlDZ1GFIlpeaw3YwDnG0bff/M3NxQpfOVJz7lCvArXz7lCvAr3xYtTORy\nHU5djpmZmbVuDwgIQEREhFwCIf9wsbOmQSCEENJAnAqau7s7BAIBqhtz1RMUE0IIIaqCU0FzdXVF\nWFgY+7qwsBB//PEHDA0NFRYYIYQQ0hCcBoW8W8wAwMzMDKNHj6b10JSIHrYmhJC6cZ7L8V1isRj3\n7t2jqa+UhB62JoSQ+nEqaHZ2djXum2lpaeGHH35QSFBEUm0PVp++9liuIyEJIUTdcS5oy5cvZ1/r\n6uqibdu2sLKyUlhgpG4vcl8Dua8BUIuNEEIAjgVt5cqV6Nmzp6JjIVIMcLZFypNXdR6TcCuDChoh\nhNc4FTQqZk3Lxc4acz27s12MRa/FVS00QgghLJkGhRDle/dh6/cHiQA0PRYhhFBBU0Pvt9jatjBG\nwq0MJNzKoAEihBDeooKmpqpbbDSknxBCqnB6sPrt27fsfI4MwyAmJgZRUVGorKxUaHCkfrR+GiGE\nVOHUQtu4cSPu3LmD8PBwHDlyBL/88gvMzMzw4MEDLF26tEFvmJycjIULF2Lu3Lnw9PRkt3/66afo\n2LEjGIaBQCBAz5494e3tze4/cOAAYmNjIRAIMGLECMyYMYPdV1hYCD8/PxQWFqKyshL+/v6wt7dv\nUFyEEELUG6eClpKSguPHj0NLSwsRERHYs2cP7O3t4eXl1aA3i4uLQ2xsLIyNjWvs69+/P9atW1fr\neVeuXMHJkycRExMDhmEwatQodO7cGW5ubgCAwMBAODg4YP78+bh+/Trmzp2Lixcv8mJF7dqG9NMA\nEUIIH3HqcjQwMICWlhbS09OhpaUFJycn6OrqwsSkYWvYODg4YNu2bTAyMmrQecePH4eHhwd0dXWh\np6eHkSNHssvWFBQU4Ny5c/j6668BAL1794auri4SEhIa9B7qqnqAiH0HC9h3sKD7Z4QQ3uLUQquo\nqMBvv/2G3377DSNHjgQA5OTkQCQSNejNWrVqJXVfeno6Zs+ejZKSEnTu3Bk+Pj6wsLAAANy9exce\nHh7ssZ07d0Z4eDiAqtajvr6+xLU7deqE5ORkDBkypEHxqStaP40QQjgWND8/P/z73/+Gubk5pkyZ\ngqysLEybNg0TJkyQWyBdunTBihUrYGhoiI0bN2LWrFk4efIkAEAoFEq0Bk1MTJCfn8/ue78L09TU\nFEKhUG6xqZN3V7tu28IYz3NLoKOjjU+7t6KiRwjRaJwKmoODg8TK1EZGRjh37pxci8aaNWvYvy9Y\nsAA9e/bEnTt30KNHDwB1Lypa277qxUj5pLYh/NXupOVRdyQhRKM16jk0Hx8fHDp0SF6xsJo1awYz\nMzNkZmaiR48eaN68OYqKitj9xcXFbHekpaUliouLJc4vKirChx9+KPX62toCmJtr3uKk/0l+We/+\nwZ90UE4wTURbW0sjP1tp+JQvn3IF+JevPEgtaDNnzkRQUBCMjY1rXT6meni9PPz1118wNjZG9+7d\nAQBlZWUoKipCy5YtAQCOjo54/Pgxe3xaWhocHR0BAPb29hCLxcjOzmaPT09Px+jRo6W+X0UFg4KC\nUrnErkrKyyvq3J/+ogB//PVEo1tp5uaGGvnZSsOnfPmUK8CvfFu0aNgAQ2mkjnL86quv2NGITk5O\n+OOPP2r8qe4ObKysrCx2kAdQtUJ2u3bt2Ot7eXnh7NmzEIvFEIlEiI2Nxfjx4wEA5ubmGD58OHu/\nLTExEeXl5eyQfj6pb7j+67fl2BmdTCteE0I0koDhcLPp1q1bcHZ25rxdmgcPHiAkJARJSUmwtbVF\nt27dsHr1amRlZSEkJARPnjxBZWUljI2N8eOPP6J9+/bsuQcPHsTp06chEAjg4eGBadOmsfvef7A6\nICAAdnZ2UuMoK6vQ2G8+7w8K+c/dLLx+Wy5xjH0HCyzy+qgpwlM4Pn2rBfiVL59yBfiVr7xaaJwK\nmjTffvstdu3aJZdAlEmTC9r7NkfcrPHgNRU0zcGnfPmUK8CvfOVV0DgNCrl//z6CgoLw5MkTiMVi\ndnteXp5cgiCKQzOJEEL4glNBW7x4MZydneHm5gZ9fX0AVYNC9u7dq9DgSOO52Flj8cSeOHvtEQDQ\n8jKEEI3FqaAZGhpi9erVNbabm5vLPSAif30cW6NbW7OmDoMQQhSKU0Hr3LkzSktLYWgo+UxETg6N\nllM37w4aodYaIUSTSC1oISEh7N+NjIwwduxY9OvXT2IKqqioqAbPuE+aDi0GSgjRZFILWkREBPr3\n78++dnR0RFFRkcSMHQ2dnJg0LWmLgVJBI4RoAqkFbcSIEVi+fHmdJwcFBck9IEIIIUQWUmcKebeY\n/e9//6ux/9dff8XQoUMVExVRiNqG67/ML4X//uvYHHGTZhAhhKg1Tgt81tYSc3Z2xpYtW+QeEFGc\n6sVA27T4Z4HV/CIRXuS+RsqTV9gZnYxjcQ+bMEJCCJFdnaMcMzMzAQBisRhZWVkSS7KYmpqivLxc\n2qlERbnYWVfdS8t9Xev+8zeeo5OtGd1XI4SonToLmru7Ozujvru7O7u9eqb9uma0J+rr13OpAEBF\njRCiVuosaHFxcWAYBr6+vti6davEPmNjY5iZ0cO66qi26bDeVT0rv76uFlqYN8MX/T6g4kYIUXmc\nJid+8eIF2rRpo4x4lIJPkxMDtU9yWv2A9cv8UuQX1f/4hXEzHZSVVwKAShc5Pk3oCvArXz7lCvAr\nX5WYbV9dUUGTlJSag1/PpdZYZqY+qvhQNp/+EwD4lS+fcgX4la/CF/gk/OFiZ41vhktfP06a2h7U\nJoSQpkIFjQCofUh/fYpei+s/iBBClETmglZayo+mMJ+42Flj9YzemOvZHc1N9es9/kXua3oYmxCi\nMjgVtJUrV9bYtnbtWgQEBMg7HqICXOyssXluP8z17A77DhZo08IIzU318f+f4JBA3Y6EEFXBafmY\np0+f1tj2008/wdfXV+4BEdXhYmctMehjc8TNGsP9qduREKIq6ixoy5YtAwA8evSI/Xs1sViMhw9p\nmiQ+qe35tepuR1Ub7UgI4R+Z76FZWFhgzZo18oyFqDgXO+taB438ei6V7qURQppcnS20devWAQCC\ng4OxYMECpQREVJupkV6NeSCrZxYZ5toW4wZ1aaLICCF8x6mFJq2Ybd++Xa7BENVX2xI01c7feE6z\n9RNCmgynQSEAEBsbi5s3b0oM17969Sq8vb0VEhhRTdXPq0mbWYRm6yeENBVOLbTt27fjypUrSEhI\ngK2tLaytrXH//n18/PHHio6PqKD6Zhb5JSYZ/vuv0301QohScSpoSUlJ2LhxI2xsbDB//nz4+Pgg\nPDwcOjqcG3hEw7jYWcNSysPXDFM1+nFndDIW77xGhY0QohScCpqenh4ASCzoqa+vzy4ASvhpnHv9\nA0CERSLsjE6mokYIUThOBY1hGDx79gxt27aFv78/EhISsGnTJrx+Xfuqx4Qfqu+n6evW/2t0+tpj\nJURECOEzTn2GM2fOxP3797FgwQLMnDkTx48fh6WlJbZs2aLo+IiKqx78sTM6uc7jXhXXv+YaIYQ0\nBqeC1rdvX/bv58+fx6tXr2Bubg5BbZP7Ed6pbqmdvvYYL/NLUV5Rc4k9C5P6JzsmhJDGkGmmEAsL\nCwgEAnz77bfyjoeoqeqZ+vcsHohhrm1r7Hfo0LwJoiKE8AmnFtr9+/cRFBSEJ0+eQCz+ZzLavLw8\nhQVG1Ff1bCHnbzxnt9HzaYQQReNU0BYvXgxnZ2e4ublBX7+q64hhGOzdu1ehwRH19Ty3pMa2ndHJ\nsDTVxzj3LlTYCCFyx6mgGRoaYvXq1TW2m5ubyz0gotmqh/EbN9PFlGFdqbARQuSG0z20zp0717pC\ndU4OPVtEalfXnI8AUPKmDDujk2nuR0KI3EhtoYWEhLB/NzIywtixY9GvXz+YmJiw26OiouDl5aXY\nCIlaql5q5kVu3c8q0r01Qoi8SG2hRUREICMjAxkZGSgqKoKjoyOKiorYbRkZGRCJ6NkiIt0X/T7g\ndFxEPLXSCCGNJ7WFNmLECCxfvrzOk4OCguQeENEc1c+nRcQ/RH6R9C8/+UUiBIYmYtW0XkqMjhCi\naaS20N4tZv/73/9q7P/1118xdOjQBr9hcnIyhg4diujoaInt9+7dg5eXFyZMmIDvvvsOhYWFEvsP\nHDiAr776CqNHj8b+/fsl9hUWFuK7777DhAkT4OXlhZSUlAbHRRTDxc4am+f2wwE/dwxzbQtpj+I/\nzS6h+2mEkEbhNCiktpaYs7Nzg6e+iouLw/79+2FsbCyxvaysDPPmzYOvry+OHj0Ke3t7+Pv7s/uv\nXLmCkydP4tixYwgPD8fJkydx+fJldn9gYCAcHBxw9OhR+Pj4YO7cuSgrK2tQbETxxg3qgm4dLKTu\nj/vfCyVGQwjRNHUWtMzMTGRmZkIsFiMrK4t9nZmZCVNTU4nZ97lwcHDAtm3bYGRkJLH9ypUr0NbW\nhqurKwBgzJgxuHjxIl69egUAOH78ODw8PKCrqws9PT2MHDkSERERAICCggKcO3cOX3/9NQCgd+/e\n0NXVRUJCQoNiI8pR1+jH8gqGWmmEEJnVWdDc3d0xaNAg3L59m/37oEGD4O7ujhEjRqBt25pTHNWl\nVatWtW6/e/cuOnbsKHGcgYEB23X4/v7OnTsjOblqMtyUlBTo6+tLXLtTp07sfqJaqu+rGejV/quX\ncCtDyRERQjRFnQ9Wx8XFgWEY+Pr6YuvWrRL7jI2NYWZmJpcghEJhjW5IU1NTCIVCdv+7jwuYmJgg\nPz+f07mx5j/lAAAgAElEQVRE9bjYWcPFzhpzNl1C2XsTGYvKKvH99qv00DUhpMHqbKHZ2tqiTZs2\n2LJlC2xtbSX+yKuYVatt5n6GYercz/VcoprcP25T6/bqh65pUVBCSENwmvqqTZva/+ORl+bNmyM1\nNVViW1FREaysrNj9RUVF7L7i4mJYWFQNLrC0tERxcXGNcz/88EOp76etLYC5uaG8wld52tpaKpnv\nnNFOiL+ZgbLyylr3h124j8GfdGjQNVU1V0XhU758yhXgX77ywKmgKZqjoyN+//139nVWVhZEIhEc\nHBzY/Y8f/7PicVpaGhwdHQEA9vb2EIvFyM7ORsuWLQEA6enpGD16tNT3q6hgUFBQcyovTWVubqiy\n+ba0aCZ1NpHi0jL4bEto0PNpqpyrIvApXz7lCvAr3xYtTOo/iAOZ1kOTt88++wzl5eVISkoCAERG\nRmLIkCHs5MdeXl44e/YsxGIxRCIRYmNjMX78eABVEyQPHz4cJ0+eBAAkJiaivLwcbm5uTZMMaZD6\nZhN5ml2CJb/8V0nREELUmYBR4s2mBw8eICQkBElJSbC1tUW3bt3YWfxTU1MREBAAbW1tmJqaYv36\n9RL36Q4ePIjTp09DIBDAw8MD06ZNY/cVFhbCz88PhYWFqKysREBAAOzs7KTGUVZWwZtvPoDqf9NL\nSs3BofOpKHkj/TGQ9i2NObXUVD1XeeNTvnzKFeBXvvJqoXEqaHl5eXj8+DFcXFxQVlaGffv2obKy\nEjNmzECzZs3kEogyUUFTTdtP3MbtdOmjU7kUNXXJVV74lC+fcgX4la9Suxw3bdqEX3/9FRUVFdi1\naxeioqJw/fp1BAYGyiUIQgDAe4wT2rc0lrqfuh8JIXXhVNCeP3+OkJAQ6OjoICoqCrt27UJYWBjS\n09MVHR/hmVXTesHKzEDq/rzCtzSbCCGkVpwKmq6uLoCqGTuaN2+OTp06AUCNB5oJkYeN3/WFgZ62\n1P3nbzzH9PXx2H7ithKjIoSoOk4FTUdHB3v27MG6deswatQoAFWz5ovFYoUGR/hr+ohu9R5zO11I\nhY0QwuJU0AIDA3H//n3Y2dlhwoQJyMrKwoYNG9ih84TIW/Wcj3VMEMOqLmyBoYmKD4wQorIaNWz/\n+fPnDZ6gWBXQKEf1kZSag53RDZtomusQf02gzp9tQ/EpV4Bf+TbJg9V5eXkSS8gsXrxYLkEQIo2L\nnTWcOlk26Jyn2SWYu/Vy/QcSQjQKp6mvYmJi8NNPP9WYM7GuCYMJkRfvMU44FvcQF/9+jsrap32s\n4a24AnO3XsbOhTRjDCF8wamg/fLLL9i/fz/s7e2ho/PPKZMnT1ZYYIS8a9ygLhg3qAsAIDA0EU+z\nS+o95624ArM2XsLeJQMVHR4hRAVw6nJs06YNevToIVHMACAkJEQhQRFSl1XTemGuZ3cYN6v/+1hF\nJYPp6+NpKRpCeIDToJCYmBhUVFRg6NChEs+eTZkyBYcOHVJogIpAg0I0S31TZlUTCIChLm3Zlp4m\n0PTP9l18yhXgV77yGhTCqcuxWbNm+PHHH7FixQp2G8MwdA+NqATvMU6cRkMyTNVD2S/zS+E9xklJ\n0RFClIVTQduyZQtWr14Ne3t7aGtXzeDAMAx8fX0VGhwhXFU/t8ZliP/tdCG2n7hNRY0QDcPpHlrr\n1q0xfPhwtGvXDra2trC1tUWbNm2wZcsWRcdHCGcudtY4td4D2lr19xzcThfSg9iEaBhOBW3QoEH4\n888/a2xft26d3AMipLH2LhkIXZ36f7WfZpfQtFmEaBBOXY6hoaHIzc2FoaEhOyiEYRgIhfXfiCek\nKexeNIDT8P7b6UIkpebAxc5aSZERQhSFU0GrXkH6XQzDUAuNqLTq6a/mbE5AWbn0J7J3RifDyswA\nG7/rq6zQCCEKwGnY/tWrV9G/f/8a2//++2/07NlTIYEpEg3b11zScq2vqAGAro4Wdi8aoKDIFIM+\nW83Fp3yVOpdj//79UVZWhri4OERGRkIsFuPp06dqWcwIP+1eNKDONdYAoKy8Et47riopIkKIvHEq\naI8ePcLQoUOxaNEi7NmzB2VlZViwYAHi4uIUHR8hcrNzoVu9g0WKS8toZhFC1BSngrZ27VosXboU\nN2/ehLW1NYyMjHD06FEcPHhQweERIl+7Fw0Al+kAdkYnU1EjRM1wKmgikQjDhw8H8M8M+8bGxtDS\natDqM4SohP1+7pyG9e+MTsaMDdRaI0RdcKpIYrEYJSWSw59LSkrw5s0bhQRFiKLtXjQA7Vsa13sc\nw1QVNnoImxDVx6mgDRs2DJ6enggJCUFubi727duHiRMn4l//+pei4yNEYVZN61XvQJFqT7NLaMAI\nISqOU0GbPn06pkyZgt9++w1ZWVk4deoUxowZgylTpig6PkIUaudCN04tNaBqwMiczQmKDYgQIjNO\nz6HV5unTp2jbtq1a3kej59A0V2Nynb4+ntNxWgJg31J3md5D3uiz1Vx8ylepz6F9++23NbZFR0dj\n4cKFcgmCEFVwwM+dUxdkJVNV/I7FPVRCVIQQrjgVtNLSmt8SvL292aVkCNEUOxe64YAft9bX+RvP\n6Zk1QlRInXM5Tp48GQKBAKmpqTXul4nFYrx9+1ahwRHSVA74uWPG+nhw6Y/fGZ0ME0NdbP++5vRw\nhBDlqbOg9e7dGwDw4sUL9OrVS2KfsbExhgwZorjICGli+/3cMXNDPCo5VLXqGUbmenanmfsJaSKc\nBoWEh4dj/PjxyohHKWhQiOZSRK6zNl5CBZeq9v8pc+Z++mw1F5/yVeqgEGnFzNfXVy5BEKLK9i4Z\nCCszA87H5xW+5TxikhAiP1K7HPfs2YOpU6dCT09P6vNmqampCguMEFVS3eJqSGtt+vp4tVyShhB1\nJbWgZWVlobKyav2o3NxczJ49W2I/wzDYu3evYqMjRMXsXTIQAPdn1srKKzF9fTznkZOEENlJLWir\nVq1i/z5lyhR8+eWXNY4RiUSKiYoQFXfAzx3eO66iuLSM0/FU1AhRPJlnClFnNChEczVFrg25XyZA\n1ehJeaHPVnPxKV+lDgohhEh3wM8d2lpcVlkDGFQVQJq9nxD5o4JGiBzsXTKwQV2KT7NLaCQkIXIm\ntcvx+PHjMDIyUtoSMcuWLUNGRgaAqgEnAoEAu3fvRrNmzQAA9+7dQ2BgILS0tGBmZob169fDzMyM\nPf/AgQOIjY2FQCDAiBEjMGPGDKnvRV2OmksVcpWlUMl6f00V8lUWPuUK8CtfhXc5hoeHs7ODhIWF\n1XpMdQGSl0OHDuHQoUMICwvDoUOH2GJWVlaGefPmwdfXF0ePHoW9vT38/f3Z865cuYKTJ0/i2LFj\nCA8Px8mTJ3H58mW5xkYIVwf83GFiqNugc6avj6cWGyGNJLWgGRgYoEWLFgCAixcv1nrMsmXLFBPV\ne65cuQJtbW24uroCAMaMGYOLFy/i1atXAKpakx4eHtDV1YWenh5GjhyJiIgIpcRGSG22f99fplYX\nzeJPiOykDtu3trbG1KlT0bp1azx69KjW4vXo0SO5BrNq1SqkpaXByMgIM2bMYOeSvHv3Ljp27Mge\n16pVKxgYGCAlJQX9+vXD3bt34eHhwe7v3LkzwsPD5RobIbJo6PB+oGoW//M3ntO8kIQ0kNSCtmnT\nJsTGxiIzMxN6enqwtbWtcYy+vr7cAunUqRM++eQTdO/eHXfv3sWUKVMQHh4OOzs7CIVCGBtLrips\namoKoVAIABAKhTAx+acP1sTEBPn5+XKLjZDGqJ6Fv6FdijujkwHIfn+NEL6RWtD09PTw1VdfsX9/\nf6aQ6u3yMnPmTPbvjo6OGDBgAI4dO8Y+4C0Q1BwW/e54ltr2E6JKDvi5Y87mBJSVVzbovOnr4+HU\nyRLeY5wUFBkhmqHO5WOqVRez/Px8ZGZmwsbGBs2bN6+1yMlL69atkZ6eDgBo3rx5jXkji4qKYGVl\nxe4vKipi9xUXF8PCwkLqtbW1BTA3N1RA1KpJW1uLN/mqeq7HfhoBAPBa+RvEZdwL2+10Iaavj0cn\nW1NsWvAZu13V85UnPuUK8C9feeBU0EpLS7F8+XKcP3+eHVI/bNgwrFmzBkZGRnIJZN++fRKtNKFQ\nCGvrqvsHjo6O+P3339l9WVlZEIlEcHBwYPc/fvyY3Z+WlgZHR0ep71VRwfBmOCzAr+G/6pLrLt8B\nAMB5EdFq6RlF+Movlu2GVJd85YFPuQL8ylepM4Vs3LgRAoEABw8exNmzZxEaGgqBQICNGzfKJQgA\nOHjwIHvf6/nz54iPj8eoUaMAAJ999hnKy8uRlJQEAIiMjMSQIUNgbm4OAPDy8sLZs2chFoshEokQ\nGxurUeu3Ec21389d5tGQNMyfEEmc5nIcP358jVGDDMNg/PjxchseHxoaigsXLkBHRwdv3rzBN998\ng5EjR7L7U1NTERAQAG1tbZiamtZ4sPrgwYM4ffo0BAIBPDw8MG3aNKnvRQ9Way51z1XWIsWHgSPq\n/tk2FJ/ylVcLjVNBmzRpEg4fPlxj+8SJE3HkyBG5BKJMVNA0lybkSkWtdprw2TYEn/JVapdjy5Yt\nsWnTJmRkZEAkEiEjIwObNm1Cy5Yt5RIEIeQfB6gbkhCZcGqhCYVCfPfdd7h79y67rXv37vjll1/Y\nkYbqhFpomksTc6UWWxVN/Gzrwqd8ldrlCFTdM7t79y4yMjJga2sLR0dHtX32iwqa5tLUXLefuI3b\n6UKZztWUwqapn600fMpX6QVNk1BB01yanuuxuIc4f+O5TOeqe2HT9M/2fXzKlwpaI1BB01x8ybUx\n98rUtbDx5bOtxqd8qaA1AhU0zcWnXAF+FTa+fbZ8ylepoxwJIarp1HqP+g+SgkZFEk3DqaDNnz8f\n27ZtU3QshBAZyDrMvxoVNqIpOBW0lJSUOmfeIIQ0PXkUtu0nbssxIkKUi1NB69atGztv4rsiIyPl\nHhAhpHEaU9iqZ/X33nFVzlERonicCtqQIUOwceNGpKamIjMzk/1z4sQJRcdHCJFRYwpbcWkZdUUS\ntcNplKOdnV3tJwsEuHfvntyDUjQa5ai5+JQr0LB8G1ucmnpUJH22mkteoxw5rYfm6uqKsLCwGtsn\nT54slyAIIYpXXZBkLWzV5zV1YSNEGk4ttLy8vFrnbCwtLYWhofqtqEotNM3Fp1yBxuWrbi02+mw1\nl9IfrE5LS8OJEyfw5s0bLFmyBOfPn8fo0aPlEoSyUUHTXHzKFZBPvupS2Oiz1VxK7XK8dOkSfH19\n4eLigqysLOjr6yMxMRE5OTn47rvv5BIIIaRpyKsr8t1rEdIUOI1y3LdvH86cOYM9e/bA3Nwcurq6\nWL9+Pa5epaG9hGiKxj7HBtBD2qRpcWqhaWlpwdbWFgDYJWMEAgG0tbUVFxkhpEk0tsX2/rnUaiPK\nwqmFVl5ejpSUFIlt6enpqKioUEhQhJCmJ48WG0CtNqI8nAaFXLlyBfPmzYOzszMeP34Me3t7/P33\n39ixYwf69eunjDjligaFaC4+5QooN9+ZG+JRKYe1OWQtkvTZai6lj3JMSUnBsWPHkJWVBRsbG4wb\nNw7dunWTSxDKRgVNc/EpV6Bp8m3M6tnva0hxo89Wc9F6aI1ABU1z8SlXoOnzlVdXIpfC1tS5Khuf\n8lVqQauoqMDu3bsRFRWFly9folWrVhg9ejRmzZqllgNDqKBpLj7lCqhOvvK8RyatuKlKrsrCp3yV\n+hzaxo0bceXKFXz55ZewsrJCbm4uYmJiUFhYiKVLl8olEEKI+nq3CDW2uNEISSIrTi00Dw8PRERE\nwNjYmN1WUlICLy8vxMbGKjRARaAWmubiU66Aaucr71abKueqCHzKV6kttJYtW0oUMwAwNjaudX5H\nQggB5PM8WzVqtREuOBU0Nzc3/PbbbxgxYgS77ezZs3B1dVVYYIQQzSDP7sj3r0HFjbxLapfjoEGD\n2L8zDIOcnBzo6+vD3NwcBQUFKC0thY2NDeLi4pQWrLxQl6Pm4lOugPrmq4gHrTWtuKnrZysLhXc5\nmpiYYPny5VJPZBgG69atk0sQhBB+kXer7f3raFpxI9xIbaFdvnwZbm5udZ7M5RhVRC00zcWnXAHN\ny5dabv/QtM+2LirxYPW3336LXbt2ySUQZaKCprn4lCugufkqau5HdSpumvrZ1kapoxxTU1MRFBSE\np0+fQiwWA6jqchQK5TP9DSGEvEsRXZLvX0udihvhhlMLbdSoURg5ciTs7e2ho1NVA6vvoUVHRys8\nSHmjFprm4lOuAL/yNTc3xFd+invuVdUKHJ8+W6W20ExNTTFz5swa2wMCAuQSBCGEcKGolltt11O1\nAkfqx6mgffTRR0hPT0enTp0ktkdFRcHZ2VkhgRFCSF0UWdzevyYVN/XAqcsxPT0dU6dOhZWVFUxM\n/mkapqamIjExUaEBKgJ1OWouPuUK8CtfrrkqazFRRRc5Pn22Su1yXLRoEYYMGYJu3bpJ3EPbu3ev\nXIIghBB5UXTLTdq1qRXX9DgVNAMDA/j7+9fYbmpqKveACCFEXpRV3KRdn4qccnEqaPb29sjPz0fz\n5s0ltqempmLw4MEKCayhxGIxVq1ahUePHqGiogI+Pj7o169fU4dFCFER7xcXZXRNSnsPKnSKwamg\nFRQU4IsvvoCzs7PEPbSrV69i/vz5CguuIXbs2AEAOHbsGJ48eYJx48bh999/r1GECSEEaJoCx+W9\nqNjJjlNBu3nzJry8vGps19PTk3tAsmAYBidPnkRISAgAoEOHDujWrRtOnz6NqVOnNm1whBC1oMzu\nybpQq052nAqal5cXZs+eXWO7hYWF3AOSxfPnz1FYWIiOHTuy2zp37oy7d+82YVSEEHVVW/FoyiJX\n/f5U1OrGqaDVVswAoH379nINRlZ5eXkAINEdamJigvT09KYKiRCiYZqyi/Ld96SiJh2ngnbjxo1a\nt2/ZsgWffvqpXANqDIFAIPG6EfMuE0JInVSxFcd3nArajBkz0KJFC7ZAFBUVQSQSwdraWqHBcWVp\naQmgKq7qQSDFxcXs9vfp6mrL7UE+dcGnfPmUK8CvfFU91zNbRtW6faRvjNzeQ9V/Bk2JU0EbPHgw\ntm7dKrHtypUreP78uUKCaqh27drBzMwMjx8/ZgtaWloaBgwY0LSBEUIIpBc6Il9aXA56v5gBwGef\nfYb4eNVoXgsEAowdOxaRkZEAgCdPniA1NRUjR45s4sgIIYQoC6cWWmZmpsRrkUiEe/fu4dmzZwoJ\nShbz58/HqlWrMG7cOFRUVGDbtm1SuxwJIYRoHk6TE9vZ2UkMuGAYBiYmJvjxxx8xahQ1pQkhhDQ9\nTi00JycniW5HXV1dWFlZQUuLU48lIYQQonCcKtKmTZtga2vL/rG2tlbpYpacnIyhQ4fWWE1bLBZj\n2bJlGDduHL7++mtcu3ZNYv/u3buxe/durFmzBmKxWJkhN4qs+f75558YPHhwjS5lVSZLrg8ePMCK\nFSuwb98+LFmyBKWl6rMkhyz5vnz5Ej/88AP27dsHb29v5OfnKztsmcj6ewwAly5dwqBBg5QVqlzI\nmu/nn3+OKVOmYMqUKUhJSVFmyDKTNddz587hwIED8PPzw++//17v+0htoc2bNw8///wzgKpRhOoi\nLi4OsbGxMDY2rrGvrvke09PT8eTJE6xbtw6nTp3CqVOnap3uS9XImi8AvH79GjY2NkqNtzFkzbWg\noADTp09Hp06dEBoaiujoaEyYMEHZ4TeYrPmWlZVh3Lhx6NOnD8LCwnDmzBl88803yg6/QRr7e3zr\n1i2lxttYjcl3zpw58PT0VGq8jSFrrkKhEImJifD394dIJGIn0KiL1GbW7du3sWzZsnr/qBoHBwds\n27YNRkZGEtur53scPXo0AMn5HgEgKSkJ9vb2AKpWF5D2MLmqkTVfoOpxDHV6+FzWXHv16sWutl5Z\nWYlmzZopN3AZyZpv27Zt0adPHwBVA7pUZUafujTm9zg0NFTt5mxtTL7x8fEIDQ1FaGgoysrKlBq3\nLGTN9erVqygvL8fBgwdx4MABtGjRot73ktpCs7W1rXUm/bdv32Lt2rX4888/VWam/Xe1atWq1u31\nzfdYUFAAKysrAIChoSEKCwsVH6wcyJqvOmpsruXl5bh16xa2bNmi0DjlpbH5btiwARkZGSo1m480\nsuZ6+/ZttGnTBhYWFmr15awxn623tzc6deqE6Oho7N+/H99++63C420MWXPNycnBy5cvsXr1avz2\n22/YvXs3FixYUOd7SW2hTZ8+XeK+ma2tLcRiMRYuXIjU1FTs378fc+fOlSW/JiFtvsfq+wvm5ubs\nvZXS0lKYmZkpP0g5qi9fTcI116CgIPj4+KjMKhGy4prv0qVL8cUXX2DTpk1KjU+e6sv1r7/+Qk5O\nDvbs2YOSkhLs3bsXFRUVTRKrPHD5bKt7G5ycnNSuq/Vd9eVqZGSErl27AqjqNUtOTq73mlIL2rBh\nwyRex8TEYPTo0TAyMsKpU6fYLg11I22+RxcXF/YGa0pKCnr16qX02BSBT/Nb1pXr/v37MXDgQHTs\n2BF//vmnskNTCGn53rhxgx3oY2Njo1aDfqSRluucOXMwe/ZszJ49G8bGxpg1axa0tbWbIkS5kpbv\nnTt32CKWmZmpVvfApZGWq7OzM7KzswFUtdZsbW3rvVa9w/ZFIhECAwMRFRWFqVOnYtGiRWr5C1Pf\nfI+dOnVChw4d8PPPP+PVq1dYsmRJk8UqD1zmtzxx4gSysrJw+PBhTJ8+ne1yVTf15frnn38iNDSU\n7dpwcnJS2y9kQP356unpITg4GB988AEePnyo8l1SdeE6T+vu3bvx+vVrhIeHY/z48UqPU17qy9fc\n3BzBwcG4fv06njx5Ah8fnyaLtbHqy9XBwQFt2rTBzz//jJcvX9bb3QjUU9DS09Ph7e2N7OxsBAcH\nY/DgwY3Noclwme9xzpw5TRSd/HHJd8yYMRgzZkwTRSg/9eXap08f/Oc//2nCCOWrvnydnJzg5OTU\nhBHKD9d5WufMmaMR/37ry7ddu3Zq3YX8Li6f7ffff9+ga0rtcoyMjMTXX38NPT09nDp1qtZipk7f\n/Pg23yOf8uVTrgC/8uVTrgC/8lVErlKnvrKzs4O2tjaGDRsGfX39Wk++evWqyn3zffDgAUJCQpCU\nlARbW1t069YNq1evBlD1EN+qVavw6NEjVFRUwNfXV627ngB+5cunXAF+5cunXAF+5avMXKUWtHHj\nxtU6y341hmHg6+uLY8eOyfzmhBBCiLxIvYc2a9asekeVzJo1S+4BEUIIIbLgNNs+IYQQoupUd4Zh\nQgghpAGooBFCCNEIVNAIIYRoBCpohBBCNAIVNEIIIRqBChohhBCNQAVNxf33v//FuHHjYGdnBy8v\nL4mZ4jds2AB3d3d8+eWXOHv2bBNGKR8//fQTPv30U4SEhHA6fteuXXB3d1fJhWa5+r//+z/281UV\nAQEBGDt2LEaPHs35s+ADrp/VF198gaSkJCVFRSQwROW9ePGCsbOzYzIzM2vsCw4OZhITE5sgKsXw\n8/NjgoODOR8fHBzM+Pn5KTAixav+fFXBX3/9xQwdOpRhGIYRi8XM4cOHmzgi1cLlsyosLJR4PWnS\npAb9ThPZUQtNjTD0DDxRsMzMTLRs2RIAoKuri4kTJzZxRKqFy79BU1NTJURCalPvemhEvZw+fRpH\njhxBs2bNIBaLMXToUEydOhUA8PTpU/z00094+/YtKioq4Obmxi65sXDhQly9ehUzZsxAWloaHj16\nhBcvXiAxMVHi+qmpqVi1ahVu376NHTt2ICoqCnfu3IGHhwe+++47rFmzBkKhEADQrFkzLF++HLa2\ntsjLy4OPjw9u3LiBtWvXIj4+HmlpaXBwcMCmTZugpVXzu9XBgwexc+dOtGvXDpMmTYKnp2etOYvF\nYqxevRoPHz5EQUEBpk6ditGjR2PTpk0ICwuDra0tfHx8MHToUCxbtgwXLlzA2LFjsXTpUonrXL16\nFVu2bEFRURFmzpyJuLg4PH36FBMmTMD06dNRWlqKOXPm4MaNGwgLC4OrqyvWrl2LqKgorFixAp6e\nnoiKisKePXtgZWUFFxcXJCUlIT8/H6tXr8bdu3eRkJCAnJwcrFy5UmISVoZhEBsbi9OnTyMrKwv2\n9vZYtWoVDA0NAVQt7Fi9bEhFRQW+/PJLdumfyZMn4/bt21i6dCmuX7+OtLQ06OvrIyoqqsbP6vXr\n19i0aRNSU1Oho6ODli1bYvny5bC0tGRjFwqFmDJlCj766KNa19sqKSnBxo0b8fDhQ+jo6EBbWxvf\nfvstPvnkEwDAmTNncOTIEejq6oJhGMydOxd9+/bl/DuwefNmJCUlwcDAAJWVlZg5cyZcXFwU9rM/\ncOAAzp8/Dz09Pejr62P58uXs2nl5eXlYuXIlXr58iZYtW2LIkCG1/g5W8/f3R3x8PLy8vDB//nz4\n+/sjNTUVmZmZSExMRM+ePeHt7c3+DB88eABdXV20adMGK1asgLGxMfbs2YOIiAh89NFHMDc3R0pK\nCt6+fYvNmzcjJiYGN2/eRFFRETZs2MB2f/7nP//Bzz//DH19fbx9+xYuLi5YtGhRnbFqpKZtIBIu\nqrs5MjIyaux7t8vx9evXjIODA5Obm8swDMPk5uYyI0aMYBiGYUpLS5kBAwYwhw4dYhiGYUQiEfPV\nV18xR48eZa81adIkZuLEiYxIJGIYhmFWrlwpNZ6uXbuy3VEPHjxgwsLCmOzsbOb06dPscf/973+Z\nadOmSZzbtWtXJigoiI3X1dWVuXDhArv/3S7HsLAwZseOHXX+bIKDgxlXV1fm2bNnDMMwzMOHDxkH\nBwfmxo0bDMMwzJIlS5hFixZJxL5w4UKp17t+/Trj6OjIXLlyhWEYhklKSmLs7e0ZoVAokcO73byT\nJk1ioqKi2NenTp1inJ2dmbS0NIZhGGbfvn1Mnz59mGvXrjEMwzC//fYb4+HhIRFT165dmQMHDjAM\nw6ps+aEAAAp9SURBVDBlZWXMpEmTmGXLljEMwzAvX75kPvroIyYuLo5hmKouLTc3N+by5cvsNQYO\nHMj4+PgwlZWVTHl5ObN69epa81uwYAHj4+PDvt6wYQMzduxYidgnT54s9efDMAwzb948iWscO3aM\n7fa9dOkS06tXLyY7O5thGIa5d+8e4+joyNy/f1/i5yftd+Dq1avM559/zh578eJFiS5lef/sf/31\nV2bQoEFMcXExwzAM8/vvvzMDBgxgxGIxwzAMM2XKFInfn9WrV9fb5fh+t3ltXY7z5s1jFixYwL4O\nDAxkvL292dfBwcFMv3792N+7wMBA5tNPP2Xz2rt3LzNnzhz2+H79+jG3bt1iGKbqZ+rm5lZnjJqK\nuhw1QPUS5jo6OjAxMcHRo0chFAphZWXFfku/dOkS8vLy4OXlBaBqVePPP/+8xmoJQ4YMgZ6eHgCw\nSzxIe88vvvgCANClSxdMmjQJlpaWyMjIwMSJEzF58mQEBQXh5s2bNc4dNmwYAMDQ0BAffPABHj9+\nXOOY4OBgPH/+nNMqtR9//DHatm0LAOjcuTNcXFzYvMeNG4fz58+jsLAQAHD8+PF6FzU1MDBA//79\nAQCOjo6oqKjA8+fP643jXR06dECnTp0AAPb29qioqEDfvn3Z10+fPpU4XiAQ4KuvvgJQ9Tl+/fXX\nOHPmDBiGQUxMDMzNzeHu7g6gqktr4MCBOH78uMQ1RowYAYFAAG1tbaxcubJGTEKhEBcuXJBY0XnC\nhAm4ffs27t69yykvoVCIP/74Q+Ianp6e7OsjR47A3d0d1tbWAKqWoXJyckJERITEdaT9DpiamiI7\nOxuxsbF4+/YtBg8eXOfvYW0a8rOPiIiAp6cnjI2NAQDDhw9HQUEBrly5guzsbFy/fh3jxo2TyLWx\nqn+GkyZNkrju+fPnUVRUxG5zdnZmF77s1q0brKysJPJ6Nw8LCwtERkYiIyMDhoaGuHDhQqPjVEfU\n5agG9PT0wDAMxGJxjX0ikYgtQHp6ejh+/DhCQ0MxatQotG/fHrNmzcKAAQOQmZkJAJgxYwZ7H+DN\nmzc17gmYmZlxjsvExETi9d69exEeHo6oqCg0b94cGRkZtS4M++49Bj09vRp5RUZGokuXLnj69Cm8\nvb3ZbjdpLCwsJF5bWlri5cuXAKqKXfv27REVFYXJkyfj77//rnfZ+nfzqv7Z1vaz53oNbW1t9j9M\noKpglZWV1Tjn3Z+9lZUVysvLkZeXh6ysLBQXF2PKlCns51VSUsIWjdrOr01GRgYEAgH7n2T1+wBV\n984cHR3rzau2a+jp6aFHjx7sdbp27SpxjqWlJfv7V03a70CPHj2wa9cuhIWF4d///jd69+4Nb29v\n9j9yLhrys8/MzMTZs2fZrnWGYWBtbY2CggJkZ2dDIBBI/H69/7smi4yMDADAli1boKenB4FAgPLy\nctja2iI3N5f92bybR/WX1Xdfv/s7efDgQYSGhmLy5MkwNzfHlClT5FJ81Q0VNDVgaWmJZs2a4enT\np+jQoYPEvidPnrDL/JSXl8PAwAD+/v5YuXIlTp48iblz5+LixYuwsbGBjo4ODh06JHF+fn6+3OK8\ndesWnJyc2P/sGloEqn3++edYvHgxJk6ciNWrV2P9+vV1Hv/q1SuJ10KhUGLpIy8vL4SFhaF169YY\nOHCgTDG9S0dHByKRiH397rfqxnj16hX7H2Zubi50dHRgZWUFGxsbWFtbS3x2lZWVKCkpadD1bW1t\nwTAMhEIhWyByc3MBADY2NjJfQywW49GjR7Czs4ONjQ17D7WaUChEly5dOF2/pKQE3bp1w44dO1Ba\nWop169Zhzpw5+OOPPwDI/2dvY2ODr776CjNnzmS3vX79Grq6uigoKADDMMjPz2dzlce/F1tbWwgE\nAqxYsYL9IgAABQUFDfpC+S6xWIxFixZh0aJFSEhIwPz589GmTRu4uLg0Ol51Ql2OakBLSwvDhg3D\nvn37JP4BJyYm4s2bN7C0tAQAZGdnY/bs2aioqIBAIMAnn3wChmHAMAwGDhwICwsLiYECJ06cwIYN\nGxocT/U13/fBBx8gJSUFpaWlAMD+J9RQhoaGEAgE2Lx5M+Li4up8xo5hGCQmJrJdgg8fPsTff/+N\n0aNHs8eMGjUK2dnZ2Lx5M9utV9f1asvtXdV5Vr/f+92H75///jWlvT5x4gQAoKysDJGRkfD09IRA\nIMCoUaOQk5OD69evs+fs2LEDBw8erDPO91laWmLo0KES3X9Hjx5Fjx492NZZfbnXdo3Dhw8jJiYG\nQFUXZnx8PLKzswFUDSK6c+eORLddXS5evIht27YBqPo9cHZ2RmVlJbtf3j/7CRMm4MyZM+yXg5KS\nEowfPx5ZWVmwtrZGv379JLp2IyMjOeXxLmNjY/bLnY+PDywtLTFkyBAcO3aMjeXhw4cS3bjvqy+v\nMWPG4PXr1wAAFxcX6OvrS/zc+EI7ICAgoKmDIPXr06cPUlJSEBQUhDNnziAyMhLPnj1DQEAAjIyM\nAFQNs37y5Al2796N06dP48yZM/jhhx/Qu3dv6OrqYuDAgQgNDcWxY8cQHR2NkpIS+Pv7Q09PDz/+\n+CMSExNx7949PH78GJ999lmtcTx79gyLFi1CTk4OEhMT0bx5c7bV2KNHD9y5cwc///wzEhMTYWBg\ngJs3byIxMRHDhg3D9OnTkZWVhdu3b6NPnz7Ytm0b/vzzTzx69Ag6OjqIj4/HuXPn2NFzOTk5SEpK\nwoULF/DgwQP2vku1Xbt24dSpU3B1dcXNmzcRGhqK06dPY+HChez9JqCqS+vZs2ewtLSssxvm77//\nxpo1a5CZmYnk5GR2hGhWVhbu3buHDz/8EK1bt0a7du2wa9cuXLx4EW/fvkV5eTmuXr0KIyMjPH36\nFDt37sTz58/x4sULWFhY4KeffmKv6eTkBF9fX/bn16ZNGwQEBCA3NxdDhw7F1q1bERoais6dO2PF\nihXQ1dWFkZER+vbti+3bt+PUqVM49f/au3sV1aEoCsAL7cTGNNpqZRtFLCTWFj6Aj5BObAQbEQQR\nSSlI0Na/ZoImlWIRFPwhjaCNgiA+QRQiBIuprtVlbnOHgcz6HuCwcziHtdkJ5OMDoVAIpVIJPp8P\nsizjdDrhcDjAtu0vu/JsNgvLsqCq6ru5aTabCAQC0DQN3W4Xt9sNy+USsVgMkUjkyzUmkwlerxcq\nlQr8fj+i0SjC4TBarRam0ylM00StVoMoinAc559nIJ1OY7FYYDAYQNM0HI9H1Ov1dx3/c+9zuRwS\niQRc14WiKDAMA4ZhQJZliKIIAJAkCbPZDL1eD/P5HKIoYrVaYbfbIZPJvO/eH9VqFaZp4nw+w3Vd\nJJNJBINBjEYjrNdrxONxpFIpSJKE7XYLVVWh6zo2mw0ajQYEQUC/38dwOMTlcoHjOLjf72i32399\nrv1+j3w+D9u20el0oOs6xuMxCoXC+x33b8IffNKvoCgKJElCOp3+6VKI6Jtw5Eiedb1eYVkWns8n\n9vs9w4zI4/hRCHnW4/FAuVyGIAgoFos/XQ4RfTOOHImIyBM4ciQiIk9goBERkScw0IiIyBMYaERE\n5AkMNCIi8gQGGhERecInh5yK783j73IAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f31646b8e10>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.semilogx(1 + np.arange(n_users), -np.sort(-user_activity), 'o')\n",
    "plt.ylabel('Number of items that this user clicked on')\n",
    "plt.xlabel('User rank by number of consumed items')\n",
    "pass"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAbwAAAE5CAYAAADrzZW0AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3XlclXX6+P/XAQ6o7KBokmiCDoKoleSUuYxKOg6MmqCo\nk2mW5ZJbk2FTJM5YLrm0DKlj6kQZEZYamOZHMytNZMZUfkSNpLjmAsJBTQ7L+f7BjzsR7sMNnnPY\nrufj4eMB9/u+z7m4gq5zv+/3ojOZTCaEEEKIJs6uvgMQQgghbEEKnhBCiGZBCp4QQohmQQqeEEKI\nZkEKnhBCiGZBCp4QQohmQQqeEEKIZkFzwSsuLmbPnj1s2bIFo9FITk6ONeMSQgghLEpTwfv55595\n5JFH+Otf/8q6desoLi7m2WefZc+ePdaOTwghhLAITQXv1Vdf5YUXXuDIkSP4+Pjg7OzM5s2b2bRp\nk5XDE0IIISxDU8ErKipi2LBhAOh0OgBcXFyws5NHgEIIIRoHTRXLaDRy7dq1SseuXbvGr7/+apWg\nhBBCCEtz0HLS0KFDGTlyJCNHjuTy5cusX7+ezz77jEcffdTa8QkhhBAWodO6W8J7771HYmIi58+f\np3379owfP54JEyYoXZxCCCFEQ6a54FUnNzcXb29vS8YjhBBCWMUdjTqZO3eupeIQQgghrEr1Gd6T\nTz7J6tWrcXFxITAwsErXpclkqnN35qlTpwgPD2fjxo2EhoYC8MMPPxAXF4ednR3u7u4sWbIEd3d3\n5ZoNGzaQkpKCTqdj+PDhTJkyRWkrKCggJiaGgoICysrKiI2NJSgoSGk/cOAAK1euxMHBgXvuuYe4\nuDgcHR3rFLsQQohGyqQiNTXVVFZWZjKZTKYxY8aYzp49W+nfmTNnTGPGjFG73Ky5c+eaevXqZUpL\nSzOZTCaT0Wg0/eEPf1C+f/PNN02zZs1Szv/qq69Mf/zjH01Go9FUVFRkGjZsmGnfvn2VXu+tt94y\nmUwm03fffWcaMGCAyWg0mkwmkyk3N9fUp08fU05OjslkMpleeOEF09KlS+sUtxBCiMZLtUtz+PDh\nyh3c888/j6+vb6V/d999NwsWLKh1gT1+/DjOzs54eXkpx7766ivs7e2Vu72oqCh2797N1atXAUhK\nSiI8PBy9Xo+joyMREREkJiYCkJ+fz86dO4mMjASgT58+6PV69u3bB0BKSgqBgYH4+fkBEBkZSXJy\nMqa6P7oUQgjRCGl6hhcXF8eVK1eqHO/Vq1et3/Ctt95i5syZlQpORkYGnTt3Vr5v164dLVq0IDMz\nEygvkre2BwQEkJGRAUBmZiZOTk60a9dOaff391faq7vWYDDIWqBCCNHMaF5pZdGiRcyZM4dvvvmm\nzm+2f/9+unTpQtu2bSsdz83NxcXFpdIxNzc3cnNzlXZXV1elzdXVlby8vDpd6+bmphwXQgjRfGia\neB4TE8OgQYO4ePEiH3/8MWvXrqVfv36MHj1a87QEk8nEv/71L+Lj46ttr24AzK13geYGyNR0rVo8\nQgghmg9NBW/QoEEAtG3blpkzZ9KzZ09efvll3nrrLY4fP67pjT777DP69etX6W6rgpeXF1lZWZWO\nGQwGWrdurbQbDAalrbCwEE9PTwC8vb0pLCyscm3Xrl2Va29tr3gdc4W6pKQUqYfVs7fXUVoqyVEj\n+VEnuTFP8qNOr7e3yOtoKnjPP/88CxYs4JNPPuHjjz/mypUrhIeHM27cOM1vlJ6ezokTJ/jmm28w\nmUzk5uby6quv0r59ex599FF27NihnHvhwgWKiooIDg4GICQkhJMnTyrtJ06cICQkBICgoCCMRiMX\nL15Uukqzs7OVQSwhISHKAJaKa93d3enYsaNqrCYT5Off0PyzNSceHq0kN2ZIftRJbsyT/Khr06bq\njVJdaHqG93//938MGDCA1NRUnnjiCfbv309cXByBgYGa32jRokVs3ryZ9957j4SEBLy9vXnxxRf5\n5z//Sb9+/SgtLSU9PR2ALVu2EBYWhoeHBwDR0dGkpqZiNBopKioiJSVFKbYeHh4MGzaM5ORkANLS\n0igpKaF///4AREREkJWVxenTp5XXjoqKkp0ehBCimdG0tFhYWBgrVqygR48ed/yGWVlZxMfH8/XX\nX9O9e3fCwsKYOHGiMvHc3t4eNze3KhPPN23axPbt29HpdISHhzN58mSl7faJ5wsXLqxUjA8ePMiK\nFStwcHCgU6dOLFq0yOzE8+LiUvmkpUI+hZon+VEnuTFP8qPOUnd4mgre999/r0xBuHr1qvL8rKmS\ngqdO/ijNk/yok9yYJ/lRZ9Muze7du/P6669z77338tBDD3HvvfeyfPlySkpKLBKEEEIIYW2aCt7b\nb7/Nf//7X15++WXWrFnDSy+9xJEjR3jrrbesHV+9WLj+O9KzLtV3GEIIISxI0yjN7777joSEBPR6\nvXIsIiKCiRMnWi2w+nTsxBWOnbjC9JHd6R3oU9/hCCGEsABNd3h6vb5SsQNwdHSscqyp2ff9ufoO\nQQghhIVoKnjOzs58+OGHlJaWAlBaWsrmzZtp2bKlVYMTQgghLEVTl+bf/vY3Hn/8cRYvXoy7uzsF\nBQX4+PiwceNGa8dXrwb28q3vEIQQQliIpoLXoUMHduzYwd69ezl//jzt27fnD3/4Q5O9w+sR0JqH\nu7eT53dCCNGEaJqHp+aDDz5gwoQJloynQbDkPLz0rEvKs8CBvXwbfRGVuULmSX7USW7Mk/yos9Q8\nPE13eGVlZezZs4eTJ09iNBqV459++mmTLHiWkp51ifitGcr3maeuyshPIYSoJ5oK3rx588jOzuZ3\nv/tdpZGZRUVFVgusKahulOe+78/VWPCa2l2hEEI0BJoK3qlTp9i2bVuVBZeb6sTz+iR3hUIIYR2a\npiV07Nix2k1Wu3XrZvGAmpLqRnnWNPJT7a5QCCHEndHcpTl79mx69+6Nm5ubcnzdunUMGTLEasE1\ndr0DfZg+srt0TwohRAOgqeCtW7eOw4cPc+XKFeztf9t59sqVK1YLrKnoHehTqyI3sJcvmaeuVjkm\nhBDizmgqeP/5z3/Yt28fTk5OlY7HxsZaJajmTO4KhRDCOjQVvICAgCrFDmD06NEWD0jU/q5QCCFE\nzTQVvNDQUBYsWMCQIUMqPcN79dVX+fTTT60WnBBCCGEpmgreihUraN26NYcOHap0PDc31ypBCSGE\nEJamqeD17NmThISEKsefeeYZiwckhBBCWMMdraVZWwkJCezZsweTycTVq1cZM2YMf/nLXwB4+OGH\n6dy5MyaTCZ1Ox/3338/s2bOVazds2EBKSgo6nY7hw4czZcoUpa2goICYmBgKCgooKysjNjaWoKAg\npf3AgQOsXLkSBwcH7rnnHuLi4nB0dFSN05JraTY1st6feZIfdZIb8yQ/6my6lqalbNmyhY0bN+Lp\n6UlOTg7Dhw+nZ8+ehISE0K9fP1577bVqr9u/fz/Jycls27YNk8nEiBEjCAgIYMCAAQDExcURHBzM\nzJkzOXToENOnT2f37t3o9Xry8vKYN28eSUlJ+Pn5ERMTw+rVq5k/f74tf3QhhBD1TNNKK5aybNky\nPD09gfLVW9zc3Dh3ruZVRJKSkggPD0ev1+Po6EhERASJiYkA5Ofns3PnTiIjIwHo06cPer2effv2\nAZCSkkJgYCB+fn4AREZGkpycjA1vbIUQQjQANi14Xbt2Vb7etWsXzs7OPPzwwwBkZ2czdepUxo8f\nT2xsLFev/jb5+vjx43Tu3Fn5PiAggIyM8vUmMzMzcXJyol27dkq7v7+/0l7dtQaDgZycHOv8kEII\nIRokmxY8gJ9++ok//elPLF26lJUrV+Li4gJAly5dWL16NZs3b8bFxYWnnnpKuSY3NxdX19/6cF1d\nXcnLy1PaKl6jgpubmzKC9PZrK6ZVyAhTIYRoXjQVvCtXrnD48GFMJhNGo5H4+Hjefvttfv3111q/\nYdeuXUlNTeWf//wn06dP5/jx4wAsXryYVq1aAfDss8+SmZnJsWPHlOuqW7zaXFtNXZbSpSmEEM2L\npkEry5cv5/r169x7772sWbOGzz77jHbt2nH27FmWLFlSpzfu1q0bAwYMICEhgWXLllVqa9myJe7u\n7pw/f54ePXrg5eWFwWBQ2gsLC5Vngd7e3hQWFla63mAwKN2nXl5eldorXsfb21s1Nnt7HR4erer0\nczV19vZ2khszJD/qJDfmSX6sT1PBO3PmDJs3bwbKdzlfv349/v7+REVFaX6j/Px80tLSeOSRR5Rj\nLVu2xGAwcOjQIZydnenevTsAxcXFGAwG2rZtC0BISAgnT55Urjtx4gQhISEABAUFYTQauXjxonJ+\ndna2MoglJCREGcBSca27uzsdO3ZUjbW01CTDg1XI0GnzJD/qJDfmSX7UWWpagqYuzYpdzo8fP46X\nlxf+/v4AVZ6dmXP9+nXi4+O5efMmAHl5eezZs4ff//73nD9/ng8//FA5NyEhAT8/P3r06AFAdHQ0\nqampGI1GioqKSElJYdy4cQB4eHgwbNgwkpOTAUhLS6OkpIT+/fsDEBERQVZWFqdPnwbKp0ZERUVV\n2cy2sUnPusTriUd4PfEI6VmX6jscIYRo8DTd4Tk4OLBu3Tr27dvHiBEjAMjIyMBoNGp+ozZt2jBo\n0CAmTZqEXq+nsLCQyMhIJkyYwIULF0hPT2fChAmUlZXh4uLCmjVrlK2I+vXrR3Z2NtHR0eh0OqKi\nopSCBuW7NsTExDB+/HjKysp45513lInlXl5erFy5knnz5uHg4ECnTp2YNWuW5rgbItkVXQghak/T\nSitnz55l1apVuLu78+KLL3L58mXmz5/P2LFjCQ8Pt0WcNtXQV1p5PfFIlT3zgjp58tfoe63+3tLt\nYp7kR53kxjzJjzqbrrRy9913s2LFCuX7u+66i4SEBLKysiwShBBCCGFtmgre+fPnqz2+cOFCZcUT\nYTuyK7oQQtSeasELDw9n8+bNuLm5MWjQIHQ6XZW5a+bmxgnrkV3RhRCi9lQL3qJFi5RVSUJDQ6vd\nHuixxx6zXmTCLNkVXQghakfzoJW77767yvH//e9/dOnSxSqB1aeGPmilPsmDdfMkP+okN+ZJftTZ\ndB5edcUO4O9//7tFghBCCCGsTdOgle+++47FixeTk5NDcXGxtWMSQgghLE5TwVu4cCGjRo0iKCgI\nJycnoHzxZbUNW4UQQoiGRlPBa926NU8//XSV44sWLbJ4QEIIIYQ1aHqGFxwczJUrV6oc379/v8UD\nEkIIIaxB9Q5vwYIFytfFxcWMGDGCe++9t9Jmql9//TUzZ860boRCCCGEBagWvEOHDvHoo48q33fq\n1KnKORULNAshhBANnWrBi46OZurUqWYvrtiEVQghhGjoVJ/h3Vrsvvjiiyrty5cvp0OHDtaJSggh\nhLAwTYNW3n///SrHJkyYUO1xIYQQoiEyOy3h8OHDABQWFpKenl5p8eiioiIMBoN1oxN1kp51SRaW\nFkKI25gteBWLQ+t0Ov7yl79UanNxcWHSpElWC0zUjeyGLoQQ1TNb8Co2eH3ssceq3S1BNDwVd3a3\nH5OCJ4Ro7jQ9w1u/fr214xBCCCGsSlPBq1g/UzR81e18LruhCyGExrU0LSUhIYE9e/ZgMpm4evUq\nY8aMUZ4N/vDDD8TFxWFnZ4e7uztLlizB3d1duXbDhg2kpKSg0+kYPnw4U6ZMUdoKCgqIiYmhoKCA\nsrIyYmNjCQoKUtoPHDjAypUrcXBw4J577iEuLq7JTpqX3dCFEKJ6mjaAtZSRI0eyceNGPD09ycnJ\nYfjw4SQmJhIYGMjQoUNZunQpoaGhvPXWW5w4cYI33ngDKF+zc8mSJWzbtg2TycSIESOIiYlhwIAB\nAMybN4/OnTszc+ZMDh06xAsvvMDu3bvR6/Xk5eUxfPhwkpKS8PPzIyYmBi8vL+bPn68ap2wAq042\nqTRP8qNOcmOe5EedTTeAVVPbWrls2TJldZaOHTvi5ubGuXPn2L9/P/b29oSGhgIQFRXF7t27uXr1\nKgBJSUmEh4ej1+txdHQkIiKCxMREAPLz89m5cyeRkZEA9OnTB71ez759+wBISUkhMDAQPz8/ACIj\nI0lOTq517EIIIRq3Oyp4jz/+eK3O79q1q/L1rl27cHZ25uGHH+b48eN07txZaWvXrh0tWrQgMzMT\noEp7QEAAGRnlQ+8zMzNxcnKiXbt2Sru/v7/SXt21BoOBnJycWsXeFKRnXeL1xCO8nniE9KxL9R2O\nEELYlOozvMGDB9d4cXVbBtXkp59+Yu7cufz666+sXr0aFxcXcnNzcXFxqXSem5sbubm5AOTm5lba\npcHV1ZW8vDylraZr27dvX6mt4nh1C2I3VTI/TwjR3KkWPFdXV1588UWgfD7eoUOHiIiIwNPTk7y8\nPLZu3crYsWNr/YZdu3YlNTWVH374gaeeeor4+HigfHL77W7tdqyu3VxbTV2Wza1LU+bnCSGaO9WC\nN3fuXB544AEA/v3vf/PWW29hZ/dbD2hYWBjTpk2rcUcFNd26dWPAgAEkJCTQvn17ZZJ7BYPBQOvW\nrQHw8vKqtIxZYWGh8izQ29ubwsLCKtdWdJ96eXlVaq94HW9vb9XY7O11eHi0qtPP1VA5ONhXe6y2\nP6e9vV2Ty40lSX7USW7Mk/xYn2rBqxgBCXDp0qVKxQ7AwcGB/Px8zW+Un59PWloajzzyiHKsZcuW\nGAwGevTowY4dO5TjFy5coKioiODgYABCQkI4efKk0n7ixAlCQkIACAoKwmg0cvHiRdq2bQtAdna2\nMoglJCREGcBSca27uzsdO3ZUjbW01NTkRks93L0dx05U7oK+y7MlL635FtA+fUFGkpkn+VEnuTFP\n8qPOpqM0W7RowerVqzl//jxFRUWcP3+eVatW0aqV9k8j169fJz4+nps3bwKQl5fHnj17ePDBB+nX\nrx+lpaWkp6cDsGXLFsLCwvDw8ADK9+ZLTU3FaDRSVFRESkoK48aNA8DDw4Nhw4aRnJwMQFpaGiUl\nJfTv3x+AiIgIsrKyOH36tPLaUVFRVQp4U1cxPy+okydBnTwZGtqBXYfPkHnqKpmnrhK/NUMGsggh\nmjRN8/BycnJ4+umnK41s9Pf355133tG8J57RaGTNmjUcOHAAvV5PYWEhYWFhzJgxAyh/Trhw4ULs\n7e1xc3OrMvF806ZNbN++HZ1OR3h4OJMnT1babp94vnDhQgIDA5X2gwcPsmLFChwcHOjUqROLFi0y\nO/G8OczDez3xCJmnrlY6FtTJk79G32v2OvkUap7kR53kxjzJjzpL3eFpnnheWlrK0aNH+eWXX7jr\nrrvo2bNnk71LkoKnTv4ozZP8qJPcmCf5UWfzief29va0atUKd3d3QkJCuHFD/sM0ZrLmphCiudG0\nlualS5eYOXMmx44d4+677+bTTz9l9OjRLFmyhPvuu8/aMQorqG7NTSi/86v4XqYsCCGaEk1dms8+\n+yy9evVi5MiRzJkzh4SEBC5evMgLL7zApk2bbBCmbTWHLs3b3T4xHWBoaAfOXL4G/FYApdvFPMmP\nOsmNeZIfdZbq0tR0h5efn6/sTlAxybtt27bNbvJ2U1bdxPRdh88oX1eszDLk951sGJUQQliOpmd4\nRqOR4uLiSseKi4v59ddfrRKUaJiqK4pCCNFYaCp4DzzwAI899hiffvopBoOBzz//nGnTpvHQQw9Z\nOz5hIzJgRQjR1Gl6hldcXMzrr79OUlISv/76Ky1btiQ6Opp58+ah1+ttEadNNcdneFD+HK/iLq5D\nG5dKXZqA0qXZHHOjlTyHUSe5MU/yo87m8/Aq5OXl4eXlZZE3b6iaa8G73a0FUAataCP5USe5MU/y\no86mg1ZudWuxe+aZZ1izZo1FAhENT+9AH5maIIRoMjQVvB9//JHVq1dz6tQpjEajcrwu++EJIYQQ\n9UFTwXv++efp1asXAwYMwMnJCSjfT+5f//qXVYMTDUdF96aDgz0Pd28nd35CiEZHU8Fr1aoVixYt\nqnK8YjcD0bTdPin92Ikrslu6EKLR0TQtISAgoNq1My9dku1kmgO13dKFEKIxUb3De/vtt5WvnZ2d\nGTNmDH379sXV9bfRMp9++inR0dHWjVAIIYSwANWCl5iYSL9+/ZTvQ0JCMBgMGAwG5VhRUZF1oxMN\nwsBevlW2EpKJ6kKIxka14A0fPpwXX3zR7MWrV6+2eECi4bl1ZwUZtCKEaKxqPfEcykdoViwi3RTJ\nxHN1MjnWPMmPOsmNeZIfdTbdADYxMZGhQ4dy+PBhADIyMoiMjOTs2bMWCUI0LulZl3g98QivJx4h\nPUsGLgkhGgdNBW/r1q3Ex8cTGhoKlD/Pi42NJS4uzqrBiYbn4PELxG/NIPPUVTJPXSV+a4YUPSFE\no6Cp4Dk6OuLv71/pWI8ePWTQSjO061BOlWMyRUEI0RhomnheVFTEtWvXcHFxUY4VFhZy8+ZNzW9U\nUlLC+++/z549e4DyHRhmz57Ngw8+CMDDDz9M586dleeD999/P7Nnz1au37BhAykpKeh0OoYPH65s\nSAtQUFBATEwMBQUFlJWVERsbS1BQkNJ+4MABVq5ciYODA/fccw9xcXE4Ojpqjl2Yl/NLIa8nHlEW\nmBZCiIZIU8EbOnQoo0ePJiIiAm9vb3Jzc/nss88YN26c5je6ePEi77//Ptu2bcPZ2ZkDBw4wffp0\ndu3ahY+PD/369eO1116r9tr9+/eTnJzMtm3bMJlMjBgxgoCAAAYMGABAXFwcwcHBzJw5k0OHDjF9\n+nR2796NXq8nLy+PefPmkZSUhJ+fHzExMaxevZr58+drjl38Zmifjhw7UXkN1es3S5QuTlmBRQjR\nUGnq0pw8eTITJkxgx44dLF26lB07dvDYY4/x+OOPa34jZ2dnZs2ahbOzMwAPPfQQTk5OHDlypMZr\nk5KSCA8PR6/X4+joSEREBImJiQDk5+ezc+dOIiMjAejTpw96vZ59+/YBkJKSQmBgIH5+fgBERkaS\nnJxMHQanCuDBkLuYPrI7QZ08cW5R9fPSv3dmyTM9IUSDpOkOT6fTMXHiRCZOnKgcy8nJqdX0BA8P\nD/785z9XOlZcXIy3tzcA2dnZTJ06lWvXrhEQEMDcuXPx9PQE4Pjx44SHhyvXBQQE8OGHHwKQmZmJ\nk5MT7dq1U9r9/f3JyMggLCyM48eP07lz50rXGgwGcnJy6NSpk6bYRWUV2wa9nnikyoT06zdLiN+a\nwd1tnHFzdpRuTiFEg6HpDu+ZZ56pcmzr1q3Mmzevzm986NAhfH196d27NwBdunRh9erVbN68GRcX\nF5566inl3Nzc3EpLmrm6upKXl6e03fpsEcDNzY3c3Nxqr3Vzc1OOiztjbrWVs5evK6M4/xr/LbHv\nHpJpDEKIeqWp4FW3cPTs2bOxt7ev05sWFRWxatUqli5dqhxbvHgxrVq1AuDZZ58lMzOTY8eOKe3m\n7iSra6upy1K6NO9cxQos1XVt3irPUFSpAMa+e0gKnxDC5sz+n+qxxx5Dp9ORlZVVqTsTwGg01mqU\n5q1iY2N54okn6NatW7XtLVu2xN3dnfPnz9OjRw+8vLwqreFZWFiodHd6e3tTWFhY6XqDwUDXrl2B\n8h3ab22veJ2KrtTq2Nvr8PBoVaeframzt7erlJshv++Es7MTyz/4j+bXOHv5OvFbM+jYzpUxg7vy\nYMhd1gi1XtyeH/EbyY15kh/rM1vw+vTpA8DZs2d54IEHKrW5uLgQFhZW6zdcsmQJPXv25JFHHsFo\nNJKbm0tOTg4uLi50794dKH+2ZzAYaNu2LVA+0f3kyZPKa5w4cYKQkBAAgoKCMBqNXLx4UTk/Oztb\nGcQSEhKiDGCpuNbd3Z2OHTuqxlhaapIlflRUt/xRtw7uylqbhutGzl6+rum1cn4pZPkH/2FoaAfG\nDu5ijXBtTpaHUie5MU/yo85SS4uZLXgzZ84Eyu+GajMFQc26desoKSlh5MiR3Lhxg/Pnz7Nz5058\nfX1JT09n8eLFACQkJODn50ePHj0AiI6OZunSpUyZMgWTyURKSgoxMTFA+WCYYcOGkZyczIwZM0hL\nS6OkpIT+/fsDEBERwZo1azh9+jR+fn5s2bKFqKgo7Ow09eYKjSoGskD50mPbvz2pufDtOnwGf193\nGdwihLCqOi0eXeHll1/m73//u6ZzT506xbBhw5TnbRUjPGfMmEFkZCRvvvkmOTk5lJWV4eLiwksv\nvVTpLmzTpk1s374dnU5HeHg4kydPVtpun3i+cOFCAgMDlfaDBw+yYsUKHBwc6NSpE4sWLTI78VwW\nj1ZXm0+h6VmXlLu+G0Ul5BnUV+ZxbuHA48MCG33Rk0/p6iQ35kl+1FnqDk9TwTMajWzYsIHvv/++\n0gCWrKws0tLSLBJIQyIFT92d/FFqufNr7BPX5X9a6iQ35kl+1Nl0t4RXX32VoqIicnJyGDVqFOHh\n4bRq1YohQ4ZYJAjRPPQO9GHRlD4MDe2geo6syymEsBZNBe/kyZPMnj2b1q1bM2rUKMaMGUN8fDzX\nrl2zdnyiCRo7uIum6QxCCGFJmgqeg0P5/5hKSkooKSkpv9DOjjNnzlgvMtGk9Q704fFhgVWOm5vM\nLoQQd0LTR2x7e3uOHTvG7373O2bMmEH//v1JT09Hr9dbOz7RhFVMXN/+7UmuFhbh6epU3yEJIZow\nTXd4s2fPpqioiFmzZlFcXMzrr7/O6dOnWbRokbXjE83A2cvXuX6zRJmQLquwCCGsQdMdnqOjI8HB\nwUD5vnRCWEp1g1T2fX+uUY/UFEI0TJru8GbOnMl///tfa8cihBBCWI2mgufk5MSBAweYMmUKixYt\n4rvvvqOsrMzasYlmoLpBKobrRunWFEJYnKaJ5z/99JOyGPOpU6fYtWsXW7dupXfv3ppXWmlMZOK5\nOmtMjlWbkN4YJ6HL5GF1khvzJD/qbDrxvGvXrphMJg4fPszmzZtJTEzk7NmzXLokn8LFnesd6IOb\nc9Wl3mQSuhDCkjQNWnnppZf48ssvuXHjBn379mXu3LkMGjSoysarQgghREOlqeC5u7vj5ubG0KFD\nGTVqlLLGwowKAAAgAElEQVQ1jxCWMrCXL5mnrlY5JoQQllKr3RJ+/PFHvvjiCzIyMujUqRNDhw7l\nvvvus2Z89UKe4amz5nOGit0VoLzYNbbndyDPYcyR3Jgn+VFnk/3wKnzxxRc88sgj/O53v6OoqIib\nN2+yZcsWEhMTOXr0qEUCEeLWPfWEEMLSNBW8N954g/T0dL744gsKCgro378/L7/8MgMHDrRyeEII\nIYRlaCp4586dIy8vj7/97W/069ePFi1aWDsuIYQQwqI0FbynnnqKGTNmWDsWIZrEczwhRMOkqeBJ\nsRO2kJ51ifitGcr3maeuNsrJ50KIhknTxHMhbEFtIWkhhLAEm205XVJSwvvvv8+ePXsAKC4uZvbs\n2Tz44IMA/PDDD8TFxWFnZ4e7uztLlizB3d1duX7Dhg2kpKSg0+kYPnw4U6ZMUdoKCgqIiYmhoKCA\nsrIyYmNjCQoKUtoPHDjAypUrcXBw4J577iEuLg5Hx6orewghhGi6bHaHd/HiRd5//33WrFlDQkIC\ns2bNYvr06Vy6dIni4mJmzJjBc889x+bNmwkKCiI2Nla5dv/+/SQnJ/PRRx/x4YcfkpyczFdffaW0\nx8XFERwczObNm5k7dy7Tp0+nuLgYgLy8PObNm8fKlStJTEzEZDKxevVqW/3Yohaqm2guk8+FEJZS\nq4KXl5dHRkYGeXl5tX4jZ2dnZs2ahbOzMwAPPfQQTk5OHDlyhP3792Nvb09oaCgAUVFR7N69m6tX\ny1feSEpKIjw8HL1ej6OjIxERESQmJgKQn5/Pzp07iYyMBKBPnz7o9Xr27dsHQEpKCoGBgfj5+QEQ\nGRlJcnIytZhvL2ykYgf0oE6e3N3GmbvbOLPv+3Oyc4IQwiI0Fbzr168zZ84c+vbtS2RkJH379mXO\nnDlcv3695ov/fx4eHvz5z3+udKy4uBgvLy+OHz9O586dlePt2rWjRYsWZGZmAlRpDwgIICOjfHBD\nZmYmTk5OtGvXTmn39/dX2qu71mAwkJOTozl2YTu9A30Y2MuXs5evc/bydTJPXZVd0IUQFqGp4C1f\nvhydTsemTZtITU1l48aN6HQ6li1bVuc3PnToEL6+voSGhpKbm1tlIWo3Nzdyc3MByM3NxdX1t6Vl\nXF1dlbvM2l7r5uamHBcNU3UDVbZ/e7IeIhFCNCWaBq38+OOPfPjhh8r3/v7+9OnTh3HjxtXpTYuK\nili1ahVLly5Vjul0uirn3drtWF271murI12ajcvZy9dJz7okUxSEEHWmqeDZ29tXOabT6ao9rkVs\nbCyTJ0+mW7duAHh5eZGVlVXpHIPBQOvWrZV2g8GgtBUWFuLp6QmAt7c3hYWFVa6t2LDWy8urUnvF\n63h7e6vGZ2+vw8OjVZ1+tqbO3t7O6rn5U9/OZJ76T5Xjqd/lMOT3naz63nfKFvlprCQ35kl+rE9T\nwWvbti3Lly9n/PjxtG7dmitXrrB582batm1b6zdcsmQJPXv2ZOjQoRiNRnJzcwkJCeHzzz9Xzrlw\n4QJFRUUEBwcDEBISwsmTv3VpnThxQtmiKCgoCKPRyMWLF5V4srOzlUEsISEhygCWimvd3d3p2LGj\naoylpSZZtVyFLVZ079bBnbvbOFfZAT3nl0L+77tTDfouT1a8Vye5MU/yo86mO56/+OKLHD58mCFD\nhtCrVy+GDBlCWloaCxYsqNWbrVu3jpKSEkaOHMmNGzc4ffo0W7ZsoX///pSUlJCeng7Ali1bCAsL\nw8PDA4Do6GhSU1MxGo0UFRWRkpKidKd6eHgwbNgwkpOTAUhLS6OkpIT+/fsDEBERQVZWFqdPn1Ze\nOyoqCjs7mXPfkP257z3VHpeJ6EKIutK8H57JZOLYsWOcP38eX19fQkJCzD5Xu92pU6cYNmyYco3J\nZEKn0zFz5kxmzJihTDy3t7fHzc2tysTzTZs2sX37dnQ6HeHh4UyePFlpu33i+cKFCwkMDFTaDx48\nyIoVK3BwcKBTp04sWrTI7MRz2Q9PnS0/hca+e6jKXV5QJ0/+Gn2vTd6/LuRTujrJjXmSH3WWusOr\n1Qawt3vjjTeYPXu2RQJpSKTgqbPlH+Xta2sCDX5tTfmfljrJjXmSH3U23QD22rVrJCcnc/LkSYxG\no3L866+/bpIFTzQMFRPRt397kquFRXi6OtV3SEKIRkzTg6y5c+fy8ccfVxkNKYQtnL18nes3Szh7\n+bpMQhdC1JmmO7xffvmFzz77rMpAj3/+859WCUqICmo7KDTkbk0hRMOk6Q6vc+fO1Y5qfPjhhy0e\nkBA1MVw31nySEELcRvUO7/Dhw8rXoaGhvPDCCzzyyCPK0lwAr776Kp9++ql1IxTN2sBevmSeulrp\n2NnL1/loz/8YO7hLPUUlhGiMVAvelClTaNOmTaUluG4tgiDrUQrr6x3oU+0k9F2HzwBI0RNCaKZa\n8Hr27ElCQoLZi5955hmLByTE7dycHeFy1Z05dh0+g7+vuzzPE0JoovoM79Zid/DgwWrPWbNmjeUj\nEuI25jaBlV0UhBBaaRq0smDBArZv3861a9esHY8QVfQO9GFoaIdq2yp2URBCiJpompbQpk0bysrK\nePHFFzGZTAwcOJDBgwcra10KYW0Vz+oqnt3dSqYpCCG00FTw1q5di5eXFyNHjuT69evs2LGD8PBw\nunbtyoYNG6wdoxBAedH7/07lVRnAItMUhBBaaOrSrNiPbuvWrcyfP59//OMfAGa32BHCGqrbRUG6\nNYUQWmi6w5s8eTKHDx/Gx8eHsLAwNmzYwH333Ver3RKEsAS1aQrbvz0p3ZpCCLM0FbywsDDs7e3x\n8PDg3nvvJSgoSIqdqDfVTVOouMuToieEUFOr7YEKCgrYs2cPX3/9NXZ2dgwcOJCIiAhrxlcvZHsg\ndQ1hC5Pqtg0CcNLbMeVPQfVa9BpCfhoqyY15kh91Nt3x/IMPPgDA0dGRli1bYjKZ2Lt3LytWrLBI\nEELURkW35u2KisuI35rB8/HfyjM9IUQVmu7wBg8eTGBgIN988w1t2rQhLCyMoUOH0qtXL1vEaHNy\nh6euoXwKVbvLu9XQ0A42X3qsoeSnIZLcmCf5UWfzDWADAgKYMWMGQUFBFnljIe6E2uCVW8l6m0KI\nW2nq0nzuueeYO3euFDvRoFQ3ReF2uw6fke5NIQSgseCNGTPG2nEIUWu9A32YPrI7Xm5OZs9L3Ps/\nG0UkhGjINHVpWlJGRgbz5s1j+vTpjBw5Ujn+8MMP07lzZ0wmEzqdjvvvv5/Zs2cr7Rs2bCAlJQWd\nTsfw4cOZMmWK0lZQUEBMTAwFBQWUlZURGxtb6W70wIEDrFy5EgcHB+655x7i4uJwdHS0zQ8srKp3\noA+9A31Iz7pE4t7/kWcoqnJOnqFI9s8TQti24O3Zs4eUlBRcXFyqtPXr14/XXnut2uv2799PcnIy\n27Ztw2QyMWLECAICAhgwYAAAcXFxBAcHM3PmTA4dOsT06dPZvXs3er2evLw85s2bR1JSEn5+fsTE\nxLB69Wrmz59v1Z9V2FZF4Xs+/ltyqyl68jxPCKGpS9NSgoODWbVqFc7OVYeUm5OUlER4eDh6vR5H\nR0ciIiJITEwEID8/n507dxIZGQlAnz590Ov17Nu3D4CUlBQCAwPx8/MDIDIykuTkZGox/VA0ImMH\nqRe0XYfP8NEe6d4UormqU8FLT08nLS2t1te1a9dOtS07O5upU6cyfvx4YmNjuXr1qtJ2/PhxOnfu\nrHwfEBBARkb5kPTMzEycnJwqvba/v7/SXt21BoOBnJycWscvGj5zWwlBedGb9cbXMpBFiGZIU8F7\n/fXXGThwIEajka1btzJp0iSmTZvGunXrLBZIly5dWL16NZs3b8bFxYWnnnpKacvNzcXV9bd5GK6u\nruTl5Sltt3eRurm5kZubW+21bm5uynHRNI0d3MVs0bv2azHxWzN44+OjNoxKCFHfNBW89PR0duzY\ngaOjI++99x5vv/02X331FTt37rRYIIsXL6ZVq1YAPPvss2RmZnLs2DGl3dzandW11dRlKV2aTdvY\nwV3wrmH05tHsXOI21r6nQgjROGkatNKiRQtatWrFuXPnyM/PZ+DAgQCV7pwsqWXLlri7u3P+/Hl6\n9OihbE9UobCwEE9PTwC8vb0pLCysdL3BYKBr165A+dZGt7ZXvI63t7fq+9vb6/DwaGWxn6cpsbe3\nazS5eSKiO8s/+I/Zc3IuXmPrt6eY9CfLzDFtTPmxNcmNeZIf69NU8IqKijh69Chbtmxh+PDhAFy/\nfh2j0TIbb3733Xe4uLjQvXt3AIqLizEYDLRt2xaAkJAQTp48qZx/4sQJQkJCAAgKCsJoNHLx4kXl\n/OzsbGUQS0hIiDKApeJad3d3s3v5lZaaZIkfFY1p+aNuHdyZPrI77+3K4tqvJarn7fruFCP7drLI\nezam/Nia5MY8yY86my4ePWPGDGbMmMGRI0eYPHkyFy5cYOjQoTzwwAMWCeLChQt8+OGHyvcJCQn4\n+fnRo0cPAKKjo0lNTcVoNFJUVERKSgrjxo0DwMPDg2HDhpGcnAxAWloaJSUl9O/fH4CIiAiysrI4\nffo0AFu2bCEqKgo7O5sOUBX1pHegD2/O7s/Q0A6o/ScvKi6T53lCNAOaFo/OysqiVatWytD+uvrp\np594++23SU9Px9fXl27durFo0SIuXLjA22+/zalTpygrK8PFxYWXXnqp0l3Ypk2b2L59OzqdjvDw\ncCZPnqy03T7xfOHChQQGBirtBw8eZMWKFTg4ONCpUycWLVpkduK5LB6trrF/Co3bmEbOxWvVtnVs\n68Irk+/sQ1xjz481SW7Mk/yos9QdnqaC1717d+bOnVtpdZOmTAqeuqbwR6k2OR3uvOg1hfxYi+TG\nPMmPOpt2ad53333VFjsZ6SgaI3OT03MuXmP6yq9knp4QTZCmghcUFMQvv/xS5fjjjz9u8YCEsLbe\ngT709FcfpXvTWEr81gwpekI0MZpGaV65coVRo0Zx3333KRO3AX7++WerBSaENc2O6mn2eR7Avu/P\n0TvQx4ZRCSGsSVPB++9//8uECROqHHdyMj+xV4iG7JXJD5gteobrlpl2I4RoGDQVvOjoaKZOnVrl\neMXkbyEaK3NF7+zl66RnXZK7PCGaCE3P8CqKncFgIDs7G5PJhMlkqvauT4jG5pXJD6guQ7b925PV\nHhdCND6aCt61a9eYM2cOffr04emnn6agoIDhw4eTnZ1t7fiEsIm2XtUv6fRLngwTF6Kp0FTwFi9e\njJeXFx999BE+Pj54eHjw5ptvsmTJEmvHJ4RNDOzlW+3xklKTLDAtRBOhqeCdPXuW2NhYevTogYND\n+WO/Ll26WGwtTSHqm7l99HIuXpOiJ0QToKngFRcXK19XTDY3mUwUFVW/WoUQjdHYwV1w0lf/JyFF\nT4jGT1PB69q1K7NmzSItLQ2j0cixY8d46aWXKq1XKURToNa1CVL0hGjsNBW8+fPnU1ZWxsSJEzl6\n9ChjxoyhsLCQ559/3trxCWFTYwd3oWNbF9X2nIvXePr1fbIKixCNkKbFoyvk5uZy7tw5fH19zW6g\n2tjJ4tHqmssCtzWtwgIwNLQDYwdXXpezueSnLiQ35kl+1Nl08eibN29y/vx5vL29CQkJ4ZtvvmHr\n1q2UlZVZJAghGppXJj9AC0d7s+fsOnyGj/b8z0YRCSHulKaCt2zZMmbNmkVxcTH//ve/efXVV4mP\nj2f58uXWjk+IevPE8G41niNFT4jGQ1PBy8zMJCkpCb1eT2JiIuvWrSM1NZW0NHmAL5qu3oE+TB/Z\nHb2Dzux5uw6fkS2FhGgENBW8Fi1aYGdnR3Z2NnZ2dvTs2RO9Xo+rq2X6VYVoqHoH+rD2r38wO5AF\nfttS6ODxCzaKTAhRW5oKXmlpKTt27GDVqlVEREQAcOnSJZmHJ5qNVyY/YHYPvQpJe36yQTRCiLrQ\nVPBiYmJ47733KCkpYeLEiVy4cIGJEyfyxz/+0drxCdFgzI7qWeOd3ulfCm0UjRCitmo1LaG5kGkJ\n6mToNLzx8VGOZueqtru20vPGrH42jKhxkN8d8yQ/6mw6LUHNc889V+trMjIyeOSRR9i6dWul4z/8\n8APR0dGMHz+eadOmUVBQUKl9w4YNPProo4wePZp33323UltBQQHTpk1j/PjxREdHk5mZWan9wIED\nREZGEh0dzYIFC2QNUHFHZkf1ZPrI7tipjGUpvFHM/HcO2DYoIUSNNG0AO3HixGqPZ2Vl1erN9uzZ\nQ0pKCi4ulbuFiouLmTFjBkuXLiU0NJS33nqL2NhY3njjDQD2799PcnIy27Ztw2QyMWLECAICAhgw\nYAAAcXFxBAcHM3PmTA4dOsT06dPZvXs3er2evLw85s2bR1JSEn5+fsTExLB69Wrmz59fq9iFuFXv\nQB+eoTvxWzOqbb9ScJM3Pj7K7KieNo5MCKFG0x3e5cuXGTVqlPJv0KBB2NnZMWnSpFq9WXBwMKtW\nrcLZ2bnS8f3792Nvb09oaCgAUVFR7N69m6tXrwKQlJREeHg4er0eR0dHIiIiSExMBCA/P5+dO3cS\nGRkJQJ8+fdDr9ezbtw+AlJQUAgMD8fPzAyAyMpLk5GSkJ1fcqd6BPrR2b6Habq7bUwhhe5oK3osv\nvlip4E2aNIn169fXegPYdu3aVXv8+PHjdO7cudJ5LVq0ULomb28PCAggI6P8k3VmZiZOTk6VXtvf\n319pr+5ag8FATk5OrWIXojrLpj2EvVrfJjD7za9tGI0QwhxNBa9fv6oP4O3t7fnf/yyzwkRubm6V\nbk43Nzdyc3OV9lvn/Lm6upKXl1ena93c3JTjQljC038OVm0rvFHME0v2ymosQjQAmp7hLViwoNL3\nRqORH3/8kbZt21osEJ2u6qfkW7sdq2vXem11pEtTWErFiixqz/OgfDWWXYfPVLvgtBDCNjQVvK+/\n/rrSXZ6zszOjRo0iKirKIkF4eXlVGQBjMBho3bq10m4wGJS2wsJCPD09AfD29qawsLDKtV27dlWu\nvbW94nXM7fZgb6/Dw6PVHfxETZe9vZ3kphpDft+JXYdPk33OYPa8XYfP8H32Fd6ZP9hGkTUc8rtj\nnuTH+jQVvMjISObMmWO1IEJCQvj888+V7y9cuEBRURHBwcFK+8mTJ5X2EydOEBISAkBQUBBGo5GL\nFy8qd5zZ2dnKIJaQkBBlAEvFte7u7nTs2FE1ntJSk8yHUSFzhdT97bHeTF3+JSWl5nsPLub9ytTX\n/o9l0x6yUWQNg/zumCf5UWfTeXjWLHYA/fv3p6SkhPT0dAC2bNlCWFgYHh4eAERHR5OamorRaKSo\nqIiUlBTGjRsHgIeHB8OGDSM5ORmAtLQ0SkpK6N+/PwARERFkZWVx+vRp5bWjoqKws7ujKYhCVGtu\n9H2azrtScFMGtAhhYzZdaeWnn37i7bffJj09HV9fX7p168aiRYuA8jl9CxcuxN7eHjc3N5YsWYK7\nu7ty7aZNm9i+fTs6nY7w8HAmT56stBUUFBATE0NBQQFlZWUsXLiQwMBApf3gwYOsWLECBwcHOnXq\nxKJFi3B0dFSNU1ZaUSefQs3z8GjF/313irXbMyjVsF2k3sGOtX8daPW4GgL53TFP8qPOUnd4srRY\nNaTgqZM/SvNuzc9He/7HrsNnNF3XHAazyO+OeZIfdVbv0kxKSiI1NdUibyJEczR2cBc2xAzStMtC\nxZ56QgjrUS14H374IQ888AAACQkJ1Z5z7tw560QlRBMyO6qnpqJ301gqRU8IK1IteC1atKBNmzYA\n7N69u9pzbp+fJ4SoXm2KXtzGNBtEJETzozotwcfHh0mTJnHXXXfx888/V1vcfv75Z6sGJ0RTMjuq\nJ+lZl8xOUAfIuXiNuI1pvDL5ARtFJkTzoHqHt3z5cv785z/j6+uLo6Mjvr6+Vf45OTnZMlYhGr3e\ngT5siBmE3sH8tJici9dkiyEhLEz1Ds/R0ZFHH31U+Xrq1KnVniOEqL21fx1Y40ayFXP1ZDNZISyj\nVtMS8vLyOH/+PO3bt8fLy8uacdUrmZagToZOm1fb/MRtTCPn4jWz5zSVuXryu2Oe5EedTVdauXHj\nBnPmzKFv375ERkbSt29f5syZw/Xr1y0ShBDN1SuTH6CFo73Zc4pLynhy6V4bRSRE06Wp4C1btgyd\nTsemTZtITU1l48aN6HQ6li1bZu34hGjy4ucNqLHolZmQoifEHdJU8H788UdWrVpFnz598Pf35/e/\n/z0rV67kxx9/tHZ8QjQL8fMG4NpKb/acMhNMWSJFT4i60lTw7O2rfvrU6XTVHhdC1M0bs/rR2r2F\n2XNMwBNS9ISoE00Fr23btixfvpxz585RVFTEuXPnWL58uUU3gBVCwLJpD9VY9EC6N4WoC02jNHNz\nc5k2bRrHjx9XjnXv3p133nlH2aS1KZFRmupkJJl5lspPTVMWADq2dWlUk9Pld8c8yY86m++WYDKZ\nOH78OOfOncPX15eQkBB0Op1FgmhopOCpkz9K8yyZHy1Fb/rI7vQO9LHI+1mb/O6YJ/lRJ9sDWZEU\nPHXyR2mepfOjZYuhxlL05HfHPMmPOpvOwxNC1I+xg7swfWR3s+es2WZ+bU4hRDkpeEI0cL0DfcwO\nZCkzwew3v7ZhREI0TlLwhGgElk17CHNPzAtvFMti00LUQFPBmzlzJqtWrbJ2LEIIM96NGWS2/UrB\nTdKzLtkoGiEaH9XdEm6VmZnJP/7xD6sGsmDBAmUHdZPJhE6nY+3atbRs2RKAH374gbi4OOzs7HB3\nd2fJkiW4u7sr12/YsIGUlBR0Oh3Dhw9nypQpSltBQQExMTEUFBRQVlZGbGwsQUFBVv15hLCGoaEd\nzA5iid+awYYaCqMQzZWmgtetWzc8PDyqHN+yZQujR4+2WDDvvfdetceLi4uZMWMGS5cuJTQ0lLfe\neovY2FjeeOMNAPbv309ycjLbtm3DZDIxYsQIAgICGDBgAABxcXEEBwczc+ZMDh06xPTp09m9ezd6\nvfmlnIRoaMYO7sLeI+coLilTPeeJJXvp6e/N7KieNoxMiIZPU5dmWFgYy5YtIysri/Pnzyv/Pv74\nY2vHB5QXNHt7e0JDQwGIiopi9+7dXL16FYCkpCTCw8PR6/U4OjoSERFBYmIiAPn5+ezcuZPIyEgA\n+vTpg16vZ9++fTaJXQhLW/vXgdjVMAX2aHYuTy370jYBCdFIaCp4MTExbNiwgZEjRzJo0CDl39Gj\nRy0azCuvvMKECROYOnUqhw4dUo4fP36czp07K9+3a9eOFi1akJmZWW17QEAAGRnlQ7UzMzNxcnKi\nXbt2Sru/v7/SLkRjtP6FmrstS8tMsu6mELfQVPBCQ0PJysqq8q93794WC8Tf35+oqCg++OADnn32\nWZ555hmysrKA8qXNXFxcKp3v5uZGbm6u0u7q+tvERFdXV/Ly8jRdK0RjVdP8vApPLNnLR3v+Z+Vo\nhGj4NBU8tRGaa9eutVggTz75JN27l/8Bh4SEMHDgQD766COlvbplzG5dJMbcMmc1XStEY9Q70Ieh\noR00nbvr8BniNqZZOSIhGjZNg1Zat27NiRMn+Pjjj/n111+ZP38+u3btsuiAldvdddddZGdnA+Dl\n5aXc7VUwGAzKwtVeXl4YDAalrbCwEE9PTwC8vb0pLCyscm3Xrl1V39veXoeHRyuL/BxNjb29neTG\nDFvn5+nRPenR1YflH/ynxnNzLl7jyWV7SX413AaRVSW/O+ZJfqxPU8H78ssvee655+jduzcXLlzA\nycmJtLQ0Ll26xLRp0ywSyPr163nyySeV73Nzc/HxKV8fMCQkhM8//1xpu3DhAkVFRQQHByvtJ0+e\nVNpPnDhBSEgIAEFBQRiNRi5evKhsZ5SdnW22WJeWmmRNOxWy3p959ZGfbh3c2RAziCeX7qWsho6L\nsjJ4NCalXnZakN8d8yQ/6my6lub69ev57LPPWLduHR4eHuj1epYsWcLXX1tuOaNNmzYpz93OnDnD\n3r17GTFiBAD9+/enpKSE9PR0oHw6RFhYmDJVIjo6mtTUVIxGI0VFRaSkpDBu3DgAPDw8GDZsGMnJ\nyQCkpaVRUlKiTFkQoqlY/8KgGndNr5Bz8ZosRyaaHU13eHZ2dvj6+gK/PQ+z9I7nU6ZMYcaMGTg4\nOPDrr78SGxurDIpxdHQkPj6ehQsXYm9vj5ubG0uWLFGu7devH9nZ2URHR6PT6YiKiqJ///5Ke2xs\nLDExMYwfP56ysjLeeecdHB0dLRa7EA3FG7P6adphAcqXI3ty6V5NIz6FaAo0bQ80btw4Xn75ZYKC\ngpg4cSLvvfce2dnZvPzyy2zevNkWcdqUbA+kTrpdzGtI+anNlARbbDHUkHLTEEl+1Nl0P7z9+/cz\nY8YMevXqxcmTJwkKCuI///kPb775Jn379rVIIA2JFDx18kdpXkPLT22Kno6a1+u8Ew0tNw2N5Eed\nzTeAzczM5KOPPuLChQu0b9+esWPH0q1bN4sE0dBIwVMnf5TmNcT8aNk5/VZ2Om0T22urIeamIZH8\nqJMdz61ICp46+aM0ryHnpy6rrlhyIeqGnJuGQPKjzqajNEtLS4mPjycsLIyQkBDCwsJYs2YNpaWl\nFglCCGF9dSleTyzZKyu1iCZD0x3ea6+9xv79+4mIiKB169ZcvnyZlJQUBg4cyAsvvGCLOG1K7vDU\nyadQ8xpDfua/c4ArBTfrdO2d3PE1htzUJ8mPOpt2aYaHh5OYmFhpTcpr164RHR1NSkqKRQJpSKTg\nqZM/SvMaU37uZGHpuhS+xpSb+iD5UWfTLs22bdtWWYDZxcVFWdpLCNH4bIgZRGv3FnW6tqKrU4jG\nRFPBGzBgADt27Kh0LDU1VdmfTgjROC2b9hAbYrSv0HK7J5bslUWpRaOh2qU5ePBg5WuTycSlS5dw\ncnLCw8OD/Px8bty4Qfv27dmzZ4/NgrUV6dJUJ90u5jX2/MRtTCPn4rU6XVtTN2djz421SX7UWapL\nUwHIZnAAABh2SURBVHVpMVdXV1588UXVC00mE6+99ppFghBCNAy3Lihd2y7LivMtOZVBCEtSvcP7\n6quvalxgWcs5jZHc4amTT6HmNbX8aF2X83bVTV5varmxNMmPugYx8fyZZ55hzZo1FgmkIZGCp07+\nKM1rqvlJz7pE/NaMOl1bccfXVHNjKZIfdTYteFlZWaxevZqcnByMRiNQ3qWZm5vL0aNHLRJIQyIF\nT538UZrX1PMz+82vKbxRXOfrpbtTXVP/3bkTNi14I0aMICIigqCgIBwcyh/7VTzD27p1q0UCaUik\n4KmTP0rzmkt+7nRKQmv3Fiyb9pCFomkamsvvTl3YtOA99thjJCQkVDn+/fff06tXL4sE0pBIwVMn\nf5TmNaf8WGoentz1lWtOvzu1ZfVRmre69957yc7Oxt/fv9LxTz/9tEkWPCFEzSoK1ZNL91J2B0vQ\nVxROa29PJISmO7zs7GwmTZpE69atcXX9rdJmZWWRltb0Jp3KHZ46+RRqXnPOj6VXXmlud37N+Xen\nJjbt0hw1ahT33nsv3bp1q/QM71//+heff/65RQJpSKTgqZM/SvMkP3c2orM6zeXOT3531Nm0S7NF\nixbExsZWOe7m5maRIIQQTUfvQB/l7swSd32mal6nud39CcvQVPCCgoLIy8vDy8ur0vGsrCyGDBli\nlcCEEI3frYXJkl2eUgBFXWjq0nzuuec4dOgQvXr1qvQM7+uvv+abb76xaoCWYDQaeeWVV/j5558p\nLS1l7ty59O3bV/V86dJUJ90u5kl+1FXk5qllX1J6J6NcaqExFUL53VFn0y7NI0eOEB0dXeW4o6Oj\nRYKwtjfffBOAjz76iFOnTjF27Fg+//zzKnesQgjr+9f8PyhfW3uLIXOvr3ewY+1fB1r1/UXDoqng\nRUdHM3Xq1CrHPT09LR6QpZlMJpKTk3n77bcB6NSpE926dWP79u1MmjSpfoMTopmzVpenFsUlZXV6\nz8Z01ygq01Twqit2AB07drRoMNZw5swZCgoK6Ny5s3IsICCA48eP12NUQojb1Wfxqw1bxyZ3opaj\nqeAdPny42uMrVqzg4YcftmhAlnblyhWASs8eXV1dyc7Orq+QhBA1uP0uqiEXQGsrLinj6df3SdGz\nAE0Fb8qUKbRp04aK8S0Gg4GioiJ8fHysGpwl6XS6St/fwSYRQggbu70AvvHxUY5m59ZTNLZXXFJW\n3yE0CZoK3pAhQ1i5cmWlY/v37+fMmdrvk2Vr3t7eQHmRrhikUlhYqByvjl5vb7FRQU2R5MY8yY86\nS+XmH9Or71mKeG6bRV6/IZLfqzunqeDdXuwA+vfvz5QpU5gwYYLFg7IkPz8/3N3dOXnypFLwTpw4\nwcCBA+s3MCGExX22YkR9hyAaME0F7/z585W+Lyoq4ocffuD06dNWCcqSdDodY8aMYcuWLdx///2c\nOnWKrKwsVqxYUd+hCSGEsCFNE88DAwMrPQMzmUy4urry0ksvMWJEw/9EdfvE8+eee44HH3ywvsMS\nQghhQ5oK3tixYyt1a+r1elq3bo2dnZ1VgxNCCCEsRVPBO336NH5+fraIRwghhLAK1Vu0GTNmKF83\ntWKXkZHBI488wtatWysdNxqNLFiwgLFjxxIZGcm3335bqX3t2rWsXbuWxYsXYzQabRmyzdQ1NwcP\nHmTIkCFVnvc2NXXJz08//cTf/vY31q9fz/z587lxo2mul1iX3Pzyyy/MmTOH9evXM3v2bPLy8mwd\ntk3U9e8K4Msvv2Tw4MG2CrVe1DU/f/zjH5k4cSITJ04kMzOzxvdRHbRy9OhRFixYUOMLvPbaazWe\n05Ds2bOHlJQUXFxcqrSZW3MzOzubU6dO8dprr/HJJ5/wySefVLu+aGNW19wAXL9+nfbt29s0Xlur\na37y8/N54okn8Pf3Z+PGjWzdupXx48fbOnyrqmtuiouLGTt2LA8++CAJCQl89tlnPP7447YO36ru\n9O/q+++/t2m8tnYn+Xn66acZOXKk5vdSvcPz9fVl5syZVf49+eSTXLp0iW3bttGhQ4fa/mz1Ljg4\nmFWrVuHs7FzpeMWam6NHjwYqr7kJkJ6eTlBQEFC+XZLa6jONWV1zA+VzNZv6ZP665ueBBx7A398f\ngLKyMlq2bGnbwG2grrnp0KGDMoDs/PnzjWK5wtq6k7+rjRs3Nvk1f+8kP3v37mXjxo1s3LiR4uLi\nGt9LteA98cQT+Pr6VvpnNBqZN28eWVlZvPvuu0yfPr2uP2O9adeuXbXHa1pzMz8/n1atWgHQqlUr\nCgoKrB+sjdU1N83FneanpKSE77//nj/96U9WjbM+3Gluli5dyrlz5xr8UoV1UdfcHD16lLvvvhtP\nT88m/WHyTn53Zs+ezeTJk/H09OTdd9+t8b1UC97QoUMrfb9t2zZGjx6Ns7Mzn3zySZMb1q+25mbF\nMwUPDw/l2cuNGzdwd3e3fZD1pKbcNHda87N69Wrmzp3baLbVsgStuXnh/7V379FUZm8cwL9KbqNS\nhzFoTTXNTIjQuNSiMRTm0oVVRibU0kUzTblUulAzhjTRbSJTp1rNQiGF0p2sEZMc5FISopJOrikZ\ninh+f1je1YlzaKbSb+zPX73v2We/z/uc3dnvu89r77VrMWvWLAQHB7/V+PpTb7m5cuUKampqwOfz\n0dTUhP3796O9vb1fYu0PfWk7XSMnenp6fRr67fUPz589ewY/Pz/Ex8dj4cKFWL16NQYPHvzKwf+/\nEDfnpqGhIQ4cOAAAuHHjBoyNjd96bP2NzUcqmaT8HDx4EBYWFvjoo4+QkZHxn7tg7I243GRlZUFD\nQwPq6upQV1f/zz/01BNxuXFzc+P2RUdHY8mSJW81rneFuPwUFBSgo6MD+vr6EAqFfXqGQGKHV1ZW\nBnd3d1RXVyMkJATTp0//F2G/23qbc3PcuHEYM2YM9uzZg4aGBnh7e/dbrG9bX+YjjY2NxYMHDxAZ\nGQlXV1coKyv3S6z9obf8ZGRk4NChQ9zQjJ6e3oDp8HrLjYyMDEJCQjB27FiUlpZi2bJl/Rbr29bX\neX737duHv//+G1FRUXB0dHzrcfaX3vKjpKSEkJAQZGZm4s6dO/D09Oy1TrEd3vHjxxEQEICxY8ci\nLi6uxwdUli1bhr179/6jk3nX9GXOzRevuAaSvuTG3t4e9vb2/RRh/+otP1OmTEF6eno/Rth/esuN\nnp4e9PT0+jHC/tPXeX7d3NwG5HdPb/n58MMPX3kIXGyH5+Pjg8GDB2PMmDEICwvrscz169df6WDv\nMjbnpngsN5Kx/IjHciMey41kbyI/YmdaeXk6sZcREVatWoWYmJh/fPD+UFJSgtDQUGRnZ0NDQwNa\nWlr45ZdfALA5N1luJGP5EY/lRjyWG8neZn7EdnjJycm9/mbXlzIMwzAM8y7o01yaDMMwDPP/ji13\nwDAMwwwIrMNjGIZhBgTW4TEMwzADAuvwGIZhmAGBdXgMwzDMgMA6PIZhGGZAYB0ewzAMMyCwDo/p\nUV5eHpydnaGpqSl2Qt/GxkYYGBjA1NQUy5cvfy3H9fPzA5/Pfy11vU7u7u4wMjJCQkJCn8oHBATA\nzMwMoaGhbziyNyc1NRW2traYNm1af4cCoHPWDXd3dzg6OsLOzg5Hjx7t75BEeHl5vVIbYd4+1uEx\nPdLX10dERATk5eWRmpqKkpKSbmXCw8MBAKamptizZ89rOe7atWvfyRWef/vtN2hqava5vK+vL6ZO\nnfoGI3rzzM3NsWHDhv4Og3Pq1Ck8fvwYUVFR4PP53VbIft3Wr1+P9evX97n8jh07XqmNMG8f6/AY\niXR0dKClpYV9+/aJ7G9pacHVq1ehq6v7Wo8nJyc3oBZIZfpOKBRCTU0NAKCiovKfXDmeebN6XQCW\nYZYuXYpVq1bBw8ODWyYqOjoaDg4OiIyM7Fb+r7/+QlhYGAYNGoS2tjY4OTlhxowZiIqKwq5duyAn\nJ4cFCxbA1dUVe/bsQUREBAwMDPD5558jPDwcKioq3N1jU1MTgoKCUFJSgiFDhmDUqFHw8fGBoqJi\nt+OmpaVh+/btaGxsxNKlS5GSkoKrV69i5cqVsLCwQFBQEJqbm9HR0QEej4eNGzdi+PDhuHnzJn76\n6Sfk5+dj165dSEhIQFlZGSwtLcVe4QcFBeHw4cPQ0dHB999/DzMzsx7LPX78GGvXrkVFRQVaWlrg\n6ekJc3NzeHt74+TJk9DU1MSmTZswadIkLFmyBLm5ufjhhx/g6uoqUk98fDz4fD6UlZUxefJkCAQC\nCIVCuLu7Y8aMGaioqMCaNWuQn5+PlJQUqKurw8vLC8nJyTh48CCMjIzA5/MRHR0NAwMDKCkp4caN\nG3j69Cm2bduGEydOIDc3F42Njdi6davInQoRITw8HKmpqaiuroapqSm8vb25haAvXbqE33//HdLS\n0ujo6MDixYthYWGBuro6eHp6IisrC7/++iuSkpJQVFQEbW3tHod66+vrERgYiOrqarS3t2P8+PHw\n9vaGgoIC+Hw+4uLi0NraChcXF1hZWcHZ2Vnk/X5+foiOjoaWlhZ8fX0xadIkeHh4cMPyly9fhp+f\nH2RkZBAeHo6ioiIcOHAAUlJSaG1thaGhIdzd3QEAISEhSEtLAwC4uLhg9OjR8Pf3R1tbG3bv3g2B\nQABZWVm0t7fDxcUFNjY2XBxCoRCrV69GeXk5ZGRksGvXLnzwwQd9as/btm1DdnY25OTk0NHRgUWL\nFsHc3LzHtsX8A8QwEjg7OxMRkY2NDW3cuJGIiFpbW2nhwoVEROTk5ERr1qzhyhcVFZGuri4VFxcT\nEVFVVRUZGRlRamoqERHx+XyytbXlyre2tnLHICKKi4sT2V6+fDmtWLGC2/bz8yN3d3ex8WZmZpKu\nri4lJSUREdHly5fp3LlzVFxczMVARBQbG0u+vr7cdmVlJY0fP55iY2OJiEgoFJKWlhYVFhZyZZyc\nnCg+Pp6IiIKDgykmJkZC5ojWrVtH06ZNo4aGBiIiSk9PpwkTJlBFRQUREc2fP5927tzJlc/KyqLg\n4GCx9cXFxZG+vj6VlJQQEdGJEyfI0NBQ5Bw0NTXp/v373D5LS0sSCATcdkhICJmamlJ9fT0RdebT\nzMyMbt26RURE+/fvJzc3N658ZmYmTZgwgU6fPk1ERE+ePKEvv/yS9uzZQ0REhYWFpKOjw+WpsrKS\nDAwMqLS0lKtj/PjxtHXrViIiamhooB07dnQ7t46ODrK3t6egoCBun4eHh8hnHRISQuvWrRObHyKi\nefPm0R9//EFERM+ePSNDQ0Oys7PjXt+wYQNVVVUREdGFCxfo7t273Gtr1qyhxMREbnvdunXdjhcY\nGEjz58+ntrY2IiJKS0sTaa9OTk60YMECam9vJyIiV1dX8vf3516X1J7T0tLoq6++4l5LSkrq9XyZ\nV8OGNJk+Wbp0KRISElBTU4P4+HjMnj27x3LR0dHQ09PDp59+CgBQVVWFpaUlDh8+DACYM2cObt26\nhYKCAgBAUlISrKyseqyrvr4eycnJcHJy4vbZ2tri/PnzaGxsFBurrKwst4rHlClTYGNjg1GjRiEz\nMxOOjo5wdnZGZGQk8vLyRN4nJSXFXamrqalh5MiRuHPnjkiZ58+fw8fHB++//z6+/fZbsTF0+eKL\nL6CkpASg87dONTU1JCYmAuhcguvYsWPo6OgAAMTExMDBwUFifaNHj8Ynn3wCoHO4uampCXV1dWLL\nUw9zw+vr63MLamppaUFZWRnjxo0DAGhra+Pu3bsi5WVlZfH1118DABQVFTFz5kzEx8dzMevr60Nb\nWxsAoKGhAUNDQ8TGxorUMWvWLACdq1T3tDL1tWvXcO3aNXz33Xfcvnnz5uH8+fN4+PChhIyIsrGx\nwfnz5wEA6enpsLe3R1FREYRCIdrb21FbWwtVVVUAgKamJkJDQ7k2kZOTg9zcXIn1Hz16FHPnzoW0\ndOfgmJmZGVasWCFSxtLSEoMGdX61amtrc22ot/Y8bNgwVFdX49SpU3j69CmmT5/OLZPDvB5sSJOR\nqOsLc9asWQgNDQWfz0dZWRkOHDjQY3mhUAhlZWWRfTweD4WFhQCAkSNHwtraGjExMZg4cSISEhLE\nLuh4//59AMD27dshIyMDKSkpPH/+HBoaGqitrcWwYcN6fF9P+7ds2YLCwkJERkZCQUEBAoGgx+HK\noUOHcv+WkZFBa2uryOv79++Hmpoabt++DWdnZ0hJSfUYQ5euzq4Lj8dDVVUVgM4v582bN+PixYsw\nMjLCkydPuCFjcV48N1lZWQDoFmNvXjxHaWnpbtsv1/fyOSgrK3Pn8ODBA9y+fRsuLi4AOtvLo0eP\nuCG8nuLuSddn/WLbUVZWBhFBKBRyHXRvrKysEBwcjLq6OiQnJ2PVqlXIyMhAUlISPv74Y5iYmHBl\nFy9eDBMTE0RFRQEAQkNDuTh68vDhQ7S0tHSLxcjISOy5ysrKcvnsrT1PnDgRe/fuRUREBPz9/WFi\nYgJ3d3fuYoT591iHx0jU9YUuLS0NV1dXBAYGwtfXl/v95mXq6uooKysT2VdfXw91dXVu28HBAW5u\nbnBwcACPxxP5wn2RhoYGpKSk4OPjg4kTJ3L7Hz16hOHDh7/SeeTl5cHMzAwKCgoAgLa2tld6f5cF\nCxbAzs4Otra2CA0N7XZ1/7KGhgaR7fr6eq4zkJGRgZ2dHaKjo1FZWQlbW9t/FFOXIUOGgIjw7Nkz\nbp+kO+G+evTokch2bW0tdw5dn+uLf0rS1tYmEkNfaGhoAADq6uq4f9fW1kJKSkqk7fSlHi0tLZw5\ncwaPHz8Gj8eDtbU1Lly4gPLycixevBhA5+dy9+5dkYue3trEyJEjIS8vj/r6epH9BQUFIu1TUmyS\n2nNTUxO0tLSwe/duNDc3Y8uWLXBzc0NycnKfz5+RjA1pMhK9OCRmb28PDw8PzJ07V2x5BwcHFBQU\noLi4GABQVVWFlJQUkaEqY2NjqKmpYeXKld2G8IiIOyaPx4OVlRViYmK4faWlpXB0dJQYb0/DeGPH\njkVeXh7a29sBABcvXuzT+16moKAAeXl5bNu2DQcOHEBOTo7EWC5evMh1GOnp6aiqqsLMmTO5Mg4O\nDsjIyEBiYqLYod0X65MUs7KyMoYOHYqioiIAnQ+T9Nbx9FYnEaG5uRmnT58G0PnQxalTpzBnzhwu\n/pycHO4ih4jg6+vLDSv2la6uLnR0dHDkyBGunpiYGNjY2PT57q6LtbU1wsLCuLs5a2tr5Obm4s6d\nO9wdtJKSEng8HrKysgAAz5494x5S6fLee+9xd2e+vr5obm7GvHnzcPz4cW7/uXPnxI52vExce+76\nv5GUlISdO3cC6Gxn+vr63HA383oM/vnnn3/u7yCYd09ZWRnc3d1RXFyMK1euYPbs2ZCWlsZnn33G\n3d0tWrQIN27cgFAoRFFREaytraGiogJdXV1s27YNCQkJOHPmDH788UdYW1uL1P/8+XOUl5eL/J4T\nFRWFQ4cO4d69e7h58yZsbGwwdepUZGZmYt++fUhMTMSVK1cQGBjY45dgTk4ONm/eDKFQiIyMDIwf\nPx4qKioAAD09PaSmpuLgwYPIzMyEoqIisrOzkZ+fj4kTJ2L16tWoqamBQCDAN998Aw8PD9y8eROl\npaVQUVHBoUOHIBAIUFRUBBUVFVy/fh3Xrl3DuXPnUFNT0+1v7gICApCSkgILCwucPXsW4eHhuHTp\nEgICAkSu7pWUlCAQCGBsbAxTU1Oxn8f58+cRFhaGe/fuobKyEpqamvDy8kJNTQ3y8/NhYmICJSUl\nqKqqYvfu3UhNTQWPx0NJSQkEAgFGjRqFtLQ0HDlyBLdv30ZzczMaGxsRGhrK1TlixAgEBARAKBQi\nPz8fQ4cORXBwMGRlZaGmpoawsDCEh4fD3NwcK1aswKBBg6CiogIdHR1s3boVCQkJOHbsGPT09ODi\n4oLm5ma4urriwYMHyM/Ph5SUFPdb38ukpKRgaWnJ5So2NhYaGhrYtGkTZGRkwOfzcezYMZSXl+PP\nP/+EsbGx2GFSHo+HgwcPwt/fH4qKihg5ciTOnj2LadOmwdDQkDuejo4OIiIicPLkSQgEAigqKiIz\nMxNPnjzB5MmTMWLECMTGxiI7Oxvy8vKwsrKCsbExKisrERISgsTERFRUVMDf3x9ycnLw9fXl2oiq\nqiqys7MREREBoVCIiooKWFpaSmzP8vLySElJwZEjRxAfH4/CwkL4+/t3Gx5m/jm24jnD9DMvLy94\nenr2+vsdwzD/DvsNj2H6QUFBAeTl5TFs2DC0tLSwzo5h3gLW4TFMP6itrcWWLVswYsQI9ug5w7wl\nbEiTYRiGGRDYU5oMwzDMgMA6PIZhGGZAYB0ewzAMMyCwDo9hGIYZEFiHxzAMwwwIrMNjGIZhBoT/\nAWwLndy7BvRhAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f3163fe9650>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.semilogx(1 + np.arange(n_items), -np.sort(-watches_per_movie), 'o')\n",
    "plt.ylabel('Number of users who watched this movie')\n",
    "plt.xlabel('Movie rank by number of watches')\n",
    "pass"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Generate co-occurrence matrix based on the user's entire watching history"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def _coord_batch(lo, hi, train_data):\n",
    "    rows = []\n",
    "    cols = []\n",
    "    for u in xrange(lo, hi):\n",
    "        for w, c in itertools.permutations(train_data[u].nonzero()[1], 2):\n",
    "            rows.append(w)\n",
    "            cols.append(c)\n",
    "    np.save(os.path.join(DATA_DIR, 'coo_%d_%d.npy' % (lo, hi)),\n",
    "            np.concatenate([np.array(rows)[:, None], np.array(cols)[:, None]], axis=1))\n",
    "    pass"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from joblib import Parallel, delayed\n",
    "\n",
    "batch_size = 5000\n",
    "\n",
    "start_idx = range(0, n_users, batch_size)\n",
    "end_idx = start_idx[1:] + [n_users]\n",
    "\n",
    "Parallel(n_jobs=8)(delayed(_coord_batch)(lo, hi, train_data) for lo, hi in zip(start_idx, end_idx))\n",
    "pass"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "User 0 to 5000 finished\n",
      "User 5000 to 10000 finished\n",
      "User 10000 to 15000 finished\n",
      "User 15000 to 20000 finished\n",
      "User 20000 to 25000 finished\n",
      "User 25000 to 30000 finished\n",
      "User 30000 to 35000 finished\n",
      "User 35000 to 40000 finished\n",
      "User 40000 to 45000 finished\n",
      "User 45000 to 50000 finished\n",
      "User 50000 to 55000 finished\n",
      "User 55000 to 60000 finished\n",
      "User 60000 to 65000 finished\n",
      "User 65000 to 70000 finished\n",
      "User 70000 to 75000 finished\n",
      "User 75000 to 80000 finished\n",
      "User 80000 to 85000 finished\n",
      "User 85000 to 90000 finished\n",
      "User 90000 to 95000 finished\n",
      "User 95000 to 100000 finished\n",
      "User 100000 to 105000 finished\n",
      "User 105000 to 110000 finished\n",
      "User 110000 to 111148 finished\n"
     ]
    }
   ],
   "source": [
    "X = sparse.csr_matrix((n_items, n_items), dtype='float32')\n",
    "\n",
    "for lo, hi in zip(start_idx, end_idx):\n",
    "    coords = np.load(os.path.join(DATA_DIR, 'coo_%d_%d.npy' % (lo, hi)))\n",
    "    \n",
    "    rows = coords[:, 0]\n",
    "    cols = coords[:, 1]\n",
    "    \n",
    "    tmp = sparse.coo_matrix((np.ones_like(rows), (rows, cols)), shape=(n_items, n_items), dtype='float32').tocsr()\n",
    "    X = X + tmp\n",
    "    \n",
    "    print(\"User %d to %d finished\" % (lo, hi))\n",
    "    sys.stdout.flush()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note: Don't forget to delete all the temporary coo_LO_HI.npy files"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "np.save(os.path.join(DATA_DIR, 'coordinate_co_binary_data.npy'), X.data)\n",
    "np.save(os.path.join(DATA_DIR, 'coordinate_co_binary_indices.npy'), X.indices)\n",
    "np.save(os.path.join(DATA_DIR, 'coordinate_co_binary_indptr.npy'), X.indptr)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.38557504805354814"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "float(X.nnz) / np.prod(X.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Or load the pre-saved co-occurrence matrix"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# or co-occurrence matrix from the entire user history\n",
    "dir_predix = DATA_DIR"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "data = np.load(os.path.join(dir_predix, 'coordinate_co_binary_data.npy'))\n",
    "indices = np.load(os.path.join(dir_predix, 'coordinate_co_binary_indices.npy'))\n",
    "indptr = np.load(os.path.join(dir_predix, 'coordinate_co_binary_indptr.npy'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X = sparse.csr_matrix((data, indices, indptr), shape=(n_items, n_items))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.38557504805354814"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "float(X.nnz) / np.prod(X.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def get_row(Y, i):\n",
    "    lo, hi = Y.indptr[i], Y.indptr[i + 1]\n",
    "    return lo, hi, Y.data[lo:hi], Y.indices[lo:hi]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "count = np.asarray(X.sum(axis=1)).ravel()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "n_pairs = X.data.sum()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Construct the SPPMI matrix"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "M = X.copy()\n",
    "\n",
    "for i in xrange(n_items):\n",
    "    lo, hi, d, idx = get_row(M, i)\n",
    "    M.data[lo:hi] = np.log(d * n_pairs / (count[i] * count[idx]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "M.data[M.data < 0] = 0\n",
    "M.eliminate_zeros()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.237642385093\n"
     ]
    }
   ],
   "source": [
    "print float(M.nnz) / np.prod(M.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now $M$ is the PPMI matrix. Depending on the number of negative examples $k$, we can obtain the shifted PPMI matrix as $\\max(M_{wc} - \\log k, 0)$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# number of negative samples\n",
    "k_ns = 1\n",
    "\n",
    "M_ns = M.copy()\n",
    "\n",
    "if k_ns > 1:\n",
    "    offset = np.log(k_ns)\n",
    "else:\n",
    "    offset = 0.\n",
    "    \n",
    "M_ns.data -= offset\n",
    "M_ns.data[M_ns.data < 0] = 0\n",
    "M_ns.eliminate_zeros()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZcAAAEdCAYAAADATSqzAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAH3hJREFUeJzt3X9QlWX+//HXkWIVOIKQxWogEbmLtGO7ts24tru4Flvt\n1NJX47CVZk47lhHlj9nVdhYbp89qn3LOLpqftdZMi5IJNpacaLaopmK18tu2C7VEnDhJmhUoHFIT\nlfv7h1/vT0fgcOLccO4bno8ZZjrXfXPfbxvxxXVd93XdLsMwDAEAYKEx0S4AADDyEC4AAMsRLgAA\nyxEuAADLES4AAMsRLgAAy4UdLg0NDcrLy1NVVVVQe3d3t1atWiWPx6N58+aprq7OPPbWW2/p+uuv\n14IFCzR//nwtWbLEusoBALZ1Vjgn1dbWaufOnUpISOh1rLS0VJJUXl4uv98vj8ejmpoaJScn67zz\nzlNpaanS0tK0Z88effLJJ9ZWDwCwpbB6Ljk5OfJ6vYqPjw9qNwxDFRUVmjt3riQpIyND2dnZqq6u\nliRNmTJFaWlpkqS///3v+vnPf25l7QAAmworXFJTU/tsb21tVWdnpzIzM822rKws1dfXB53X3d2t\nkydPaty4cRGUCgBwiogm9Nva2iRJbrfbbHO73Tp48GDQea+88opmz54dya0AAA5iydNiLpcr6POZ\n25XV1dVp1qxZVtwKAOAAEYVLSkqKJCkQCJhtXV1dZrskdXR0KCkpqVcAAQBGrrCeFutPenq6EhMT\n1dLSouTkZElSc3OzcnNzzXOSkpK0bNmysK534sRJOW2P5pgYl06edFbRTqvZafVK1DwcnFav5Mya\nzz47ZlDfF1G4uFwuFRQUqLKyUjNmzJDf71djY6PWr18/qOsZhtTRcSSSkoZdUlIcNQ8xp9UrUfNw\ncFq9kjNrnjjRPfBJfQhrWKypqUnFxcXy+XwqKytTSUmJeayoqEiGYcjj8WjFihXyer1Bw2IAgNHH\nZaeXhR0/ftJxqe7E30ScVrPT6pWoeTg4rV7JmTUPac8FAIBvIqI5F6v96ZEd+qC1q89j/+dnF+vS\n708f5ooAAINhq3A52HVCB45/u89jbe0H+2wHANiPrcKlP4bRo7Yv2tXaurfP45Mnn68xYxjhAwC7\ncES4HO1qU8U/2rWz/h+9jn3V1a4Nv7lWaWnpUagMANAXR4SLJI11pyhu/LnRLgMAEAbHhEt/DKNH\nn366v89jDJcBQHQ4Ply++vKQ/vDE/9VYtz+4neEyAIgax4eLxJAZANgNY0YAAMuNiJ5LX0LNxUjM\nxwDAUBqx4dLfXIzEfAwADLURGy4SczEAEC2MCwEALDeiey79YT4GAIbWqAwX5mMAYGiNynCRmI8B\ngKHE2A8AwHKjtufSH/YqA4DIES5nYK8yAIgc4dIH5mMAIDKM8QAALEfPJUz9zcV0do5TIHCU+RgA\n+BrCJUysjQGA8IUdLg0NDVq2bJmWLFmi/Px8s727u1urV6/WRx99pJMnT2rp0qWaNWuWefyFF17Q\n/v371dTUpJ/+9Ke6+uqrrf0TDCPmYgAgPGGFS21trXbu3KmEhIRex0pLSyVJ5eXl8vv98ng8qqmp\nUXJystrb2/XWW2+ppKREx44dU1tbm7XVAwBsKaxJgpycHHm9XsXHxwe1G4ahiooKzZ07V5KUkZGh\n7OxsVVdXS5Jef/11nThxQo8//rgee+wxTZw40eLyAQB2FFa4pKam9tne2tqqzs5OZWZmmm1ZWVmq\nr6+XJH3++ec6cOCAFi5cqClTpmjz5s0WlGw/pyf7W1v39vnV09MT7RIBYFhFNKF/epjL7XabbW63\nWz6fT5IUHx+v73znO5KkadOm6W9/+1skt7MtJvsBIJglT4u5XK6gz4ZhSJIuueQSbdu2TdKpXszk\nyZOtuJ0tMdkPAP8ronBJSUmRJAUCASUnJ0uSurq6zPacnBydf/75evjhh3XgwAHdddddEZYLAHCC\niMIlPT1diYmJamlpMcOlublZubm55jnFxcXhX3CMa+BzHMYwetTVdVCdneN6HUtLS4vKwsuYmDFK\nSoob9vsOltPqlah5ODitXsmZNQ9WROHicrlUUFCgyspKzZgxQ36/X42NjVq/fv3gLthjRFKOLX31\n5SHd+z91GutOCW6P4lxMUlKcOjqODPt9B8tp9UrUPBycVq/kzJonTnQPfFIfwvq1uampScXFxfL5\nfCorK1NJSYl5rKioSIZhyOPxaMWKFfJ6veawGE45PR/z9a8zwwYARpKwei5Tp041F0ueKTY2VmvX\nrrW0KACAs7HTIgDAcmxcGSWh3ngp8dZLAM5GuEQJCy8BjGSESxSx8BLASMW4CwDAcvRcbCjUfAxz\nMQCcgHCxof7mY5iLAeAUhItNMR8DwMkYXwEAWI5wAQBYjmExB2HhJQCnIFwchIWXAJyCcHEYJvoB\nOAFjKAAAyxEuAADLMSw2QjDZD8BOCJcRgsl+AHZCuIwgTPYDsAvGSQAAliNcAACWY1hsFDhzsr+z\nc5wCgaOSmOgHMDQIl1GALfwBDDfCZZRgsh/AcGI8BABgubDDpaGhQXl5eaqqqgpq7+7u1qpVq+Tx\neDRv3jzV1dUFHb/66qu1YMECLViwQO+//741VQMAbC2sYbHa2lrt3LlTCQkJvY6VlpZKksrLy+X3\n++XxeFRTU6Pk5GRJ0uLFi5Wfn29hybAKq/oBDJWwwiUnJ0dz5szR/Pnzg9oNw1BFRYU2btwoScrI\nyFB2draqq6u1cOFCSdLLL7+sQ4cOSZJuvvlmnX322RaWj0iwqh/AUAnr19LU1NQ+21tbW9XZ2anM\nzEyzLSsrS/X19ebnu+++W7feeqsmTJigLVu2RFgurHZ6ov/Mr7HulGiXBsDBIhrzaGtrkyS53W6z\nze126+DBg+bnCy+8UJI0ffp0vfvuu5HcDgDgEJYMqLtcrqDPhmFIkv7973+bgbJ//35NmjTJitsB\nAGwuonUuKSmnhk4CgYA5gd/V1WW2JyUlacOGDXrzzTfl9/u1dOnSCMvFcGGyH0AkIgqX9PR0JSYm\nqqWlxQyX5uZm5ebmmscffPDB8C84xjXwORgWA032b7u/QFOmTBmWWmJixigpKW5Y7mUVah56TqtX\ncmbNgxVRuLhcLhUUFKiyslIzZsyQ3+9XY2Oj1q9fP7gL9hiRlAOLhVrVHwgcVUfHkWGpIykpbtju\nZRVqHnpOq1dyZs0TJ7oHPqkPYY1rNDU1qbi4WD6fT2VlZSopKTGPFRUVyTAMeTwerVixQl6v1xwW\nAwCMTmH1XKZOnWouljxTbGys1q5da2lRAABnY+NKfGOhJvuZ6AcgES4YBLbwBzAQwgWDwhb+AEJh\n/AIAYDl6LrAMCy8BnEa4wDLssgzgNMIFlmIuBoDEnAsAYAgQLgAAyzEshmHBZD8wuhAuGBZM9gOj\nC+GCYcNkPzB6MA4BALAc4QIAsBzDYoi6gXZZBuA8hAuibqBdlpOTvxudwgAMGuECW2CyHxhZmHMB\nAFiOcAEAWI5hMdjW6Yn+8ePHKRA42us4q/oB+yJcYFus6geci3CBrTHRDzgTYwoAAMvRc4EjDbTw\nkrkYILoIFzjSQAsvmYsBoivsX+8aGhqUl5enqqqqoPbu7m6tWrVKHo9H8+bNU11dXa/vfeWVVzRn\nzpzIqwW+5vR8zNe/xrpTol0WAIXZc6mtrdXOnTuVkJDQ61hpaakkqby8XH6/Xx6PRzU1NUpOTpYk\nHT58WO+++66FJQMA7C6snktOTo68Xq/i4+OD2g3DUEVFhebOnStJysjIUHZ2tqqrq81ztm7dqoUL\nF1pXMQDA9sLquaSmpvbZ3traqs7OTmVmZpptWVlZqq+vlyT961//0vnnn68JEybIMAwLygVC43XK\ngD1ENKHf1tYmSXK73Wab2+2Wz+eTJO3evVsul0uPPPKIvvzySz366KNatGiRYmJiIrkt0C8WXgL2\nYMnTYi6XK+jz6V7K4sWLzbYdO3bo17/+tRW3A0Ji4SUQfRGFS0rKqSdzAoGAOYHf1dVltp+2efNm\nHT58WE8//bR+9atf9X/BMa7+jwEWGD9+nJKS4ob8PjExY4blPlZyWs1Oq1dyZs2DFVG4pKenKzEx\nUS0tLWa4NDc3Kzc3N+i8xYsXB/Vi+tXDvAyGViBwVB0dR4b8PklJccNyHys5rWan1Ss5s+aJE90D\nn9SHiGY2XS6XCgoKVFlZKUny+/1qbGzUtddeG8llAQAOF1bPpampSRs3bpTP51NZWZneeecdrVmz\nRpJUVFSk1atXy+Px6OTJk/J6vb2GxQA74EkyYPiEFS5Tp041F0ueKTY2VmvXrrW0KGAo8CQZMHzY\nWwyjCk+SAcODMQAAgOUIFwCA5QgXAIDlmHMBxMvHAKsRLoB4+RhgNcIF+P94kgywDuEChMDCS2Bw\nCBcgBBZeAoNDuAADYLgM+ObozwMALEfPBRikUPMxnZ3j5HanMB+DUYtwAQaJ+Rigf4QLEAHmY4C+\n0WcHAFiOcAEAWI5wAQBYjjkXYAiwESZGO8IFGAJshInRjnABhghPkmE0o28OALAc4QIAsBzhAgCw\nHOECALBc2OHS0NCgvLw8VVVVBbV3d3dr1apV8ng8mjdvnurq6sxjBw4c0D333KO//OUvuvvuu3Xw\n4EHrKgcA2FZY4VJbW6stW7YoISGh17HS0lJJUnl5uR566CEtW7bMDJHjx4/L4/Hotttu06WXXqrn\nnnvOwtIBAHYVVrjk5OTI6/UqPj4+qN0wDFVUVGju3LmSpIyMDGVnZ6u6ulqSlJaWppkzZ0qS9u/f\nrylTplhZOwDApsJa55Kamtpne2trqzo7O5WZmWm2ZWVlqb6+Pui8Bx54QPv27dPll18eQamA84Va\nuS+xeh8jR0SLKNva2iRJbrfbbHO73fL5fEHn/fa3v9VLL72kBx98UKtWrYrkloCj8Q4YjBaWrNB3\nuVxBnw3DkCS9/fbbmjx5siZNmqRJkyZp//7+f2MDRgtW7mM0iChcUlJSJEmBQEDJycmSpK6uLrM9\nNjZWGzZs0AUXXKAPP/xQt99+e+gLjnGFPg6McOPHj1NSUlxU7h0TMyZq9x4Mp9UrObPmwYooXNLT\n05WYmKiWlhYzXJqbm5WbmytJmj59uqZPnx7+BXuMSMoBHM0wevTBBz4FAkf7PD7U8zFJSXHq6Dgy\nZNe3mtPqlZxZ88SJ7oFP6kNE4eJyuVRQUKDKykrNmDFDfr9fjY2NWr9+fSSXBUYl5mMwkoQVLk1N\nTdq4caN8Pp/Kysr0zjvvaM2aNZKkoqIirV69Wh6PRydPnpTX6zWHxQB8M8zHYKQIK1ymTp1qLpY8\nU2xsrNauXWtpUQAAZ+N9LoAD8GZLOA3hAjgAb7aE0xAugEMwHwMnoS8NALAc4QIAsBzhAgCwHOEC\nALAc4QIAsBzhAgCwHI8iAw7Gy8dgV4QL4GBsdgm7IlwAh2NxJeyIcAFGKPYjQzQRLsAIxX5kiCbC\nBRjB+hoy669H09k5ToHAUXo1sAThAowyPASA4UC4AKMQDwFgqNH3BQBYjnABAFiOcAEAWI5wAQBY\njnABAFiOcAEAWI5HkQGEpaenR/v2fdLnMRZe4kzfKFwaGhq0bNkyLVmyRPn5+WZ7d3e3Vq9erY8+\n+kgnT57U0qVLNWvWLElSU1OTtm3bpgsuuEBNTU267777FBcXZ+2fAsCQ27fvE931389prDslqJ2F\nl+hL2OFSW1urnTt3KiEhodex0tJSSVJ5ebn8fr88Ho9qamqUnJysjo4OLVq0SBdeeKG2bt2qqqoq\n3Xjjjdb9CQAMGxZfIlxhh0tOTo7mzJmj+fPnB7UbhqGKigpt3LhRkpSRkaHs7GxVV1dr4cKFuuyy\ny8xze3p6NG7cOItKB2C1UDsph3opGXCmsMMlNTW1z/bW1lZ1dnYqMzPTbMvKylJ9fX3QeSdOnNC7\n776r9evXD7JUAEMt1L5jHQealZSaNfxFwZEintBva2uTJLndbrPN7XbL5/MFnffHP/5RS5cuVWxs\nbKS3BDCE+hv6OtrVHoVq4FSWPd7hcrmCPhuGYf73li1bNHv2bGVmZmrXrl1W3RIAYFMR91xSUk49\nORIIBJScnCxJ6urqMtt37dqlrVu3msNm06dP18yZM/u+2BhX3+0AbMswetTVdVCdnX3Pp6alpfX5\nmHJPT49aW1u/0fd8XUzMGCUlOevJUyfWPFgRh0t6eroSExPV0tJihktzc7Nyc3MlSTNnztQbb7wR\n3sV6jIHPAWArX315SPf+T12vR5Sl0I8pt7bujejR5qSkOHV0HIms+GHmxJonTnQPfFIfIg4Xl8ul\ngoICVVZWasaMGfL7/WpsbGTiHhhFBvuIMo82j1xhz7k0NTWpuLhYPp9PZWVlKikpMY8VFRXJMAx5\nPB6tWLFCXq/XHBYDAIw+Yfdcpk6dai6WPFNsbKzWrl1rWVEAAGdjMyAAgOUIFwCA5QgXAIDlCBcA\ngOUIFwCA5XhZGIAhwy7LoxfhAmDIsMvy6EW4ABhS33SX5VC9HYlXKjsF4QLAVkL1dnilsnMQLgBs\nhz3HnI++JQDAcoQLAMByhAsAwHKECwDAcoQLAMByhAsAwHKECwDAcoQLAMByhAsAwHKECwDAcmz/\nAmBE6+np0b59n/R7nI0whwbhAmBE27fvE931389prDul1zE2whw6hAuAEY+NMIcffUEAgOW+Ubg0\nNDQoLy9PVVVVQe3d3d1atWqVPB6P5s2bp7q6uqDju3bt0hVXXKH9+3mtKQCMBmEPi9XW1mrnzp1K\nSEjoday0tFSSVF5eLr/fL4/Ho5qaGiUnJ0uSDh8+rEmTJllUMgDA7sIOl5ycHM2ZM0fz588PajcM\nQxUVFdq4caMkKSMjQ9nZ2aqurtbChQslSVdccYW2bdtmXdUAYIFQr1QerqfIRurTbGGHS2pqap/t\nra2t6uzsVGZmptmWlZWl+vr6yKsDgCHU3yuVh/MpspH6NFvET4u1tbVJktxut9nmdrvl8/kivTQA\nDDk7PElmhxqsZtmjyC6XK+izYRhWXRoAJAUPY3V2jlMgcNQ81tPTI0m9hpD6G/YK9z59cepQ1XCK\nOFxSUk515QKBgDmB39XVZbZL0jPPPKNPP/1UTz75pBYtWqRzzjmn74uNcfXdDgDqfxhLkjoONGts\n/IRew0sdB5qVlJpl2X2+6mrXtvsLNGXKlG90TUmKiRmjpKS4oLbOznEhv2f8+HG9vscJIg6X9PR0\nJSYmqqWlxQyX5uZm5ebmmufccMMNuuGGGwa+WA+9HQCh9TeEdLSrvc9jR7vaLb2PJAUCR9XRceQb\nXzMpKa7X932992XlvawycaJ74JP6EHG/zuVyqaCgQJWVlZIkv9+vxsZGXXvttZFeGgDgUGGHS1NT\nk4qLi+Xz+VRWVqaSkhLzWFFRkQzDkMfj0YoVK+T1eoOGxQAAo0vYw2JTp041F0ueKTY2VmvXrrWs\nKACAs/G4AwDAcoQLAMBybLkPAEPs9BYvZ67NkQa3DscJCBcAGGKhtngZzDocJyBcAGAYhFqfMxIx\n5wIAsBzhAgCwHOECALAccy4AYFN2eJnZYBEuAGBTdniZ2WARLgBgY059kZh9+1QAAMciXAAAliNc\nAACWI1wAAJYjXAAAliNcAACWI1wAAJYjXAAAlmMRJQBY5PRLwc40Ul8IFgrhAgAW6e+lYCP1hWCh\nEC4AYKG+tmsZqS8EC4U5FwCA5QgXAIDlwg6XhoYG5eXlqaqqKqi9u7tbq1atksfj0bx581RXVxd0\nfPPmzdq8ebP+67/+S93d3dZUDQCwtbDCpba2Vlu2bFFCQkKvY6WlpZKk8vJyPfTQQ1q2bJkOHjwo\nSfL5fPL7/Vq8eLGys7P117/+1cLSAQB2FVa45OTkyOv1Kj4+PqjdMAxVVFRo7ty5kqSMjAxlZ2er\nurpakrRnzx5NmzZNkjRt2jS9/fbbVtYOALCpsMIlNTW1z/bW1lZ1dnYqMzPTbMvKylJ9fb0kqaOj\nQ3FxcZKkuLg4dXZ2RlovAMABIprQb2trkyS53W6zze12m8NiSUlJOnLkiCTpyJEjSkxMjOR2AACH\nsORpMZfLFfTZMAxJ0qWXXqr3339fkvT+++/rsssus+J2AACbcxmnkyAM8+fP19y5c5Wfny9J+vjj\nj3XVVVeprq5OycnJkqT7779fhw4d0vr16yWdelrsxIkTOnTokH7zm98oNjZ2CP4YAAA7iWiFfnp6\nuhITE9XS0mKGS3Nzs3Jzc81zFi9eHFGBAADniWhYzOVyqaCgQJWVlZIkv9+vxsZGXXvttZYUBwBw\nprCGxZqamrRx40bt2bNHkydPVnZ2ttasWSPp1CLK1atX66OPPtLJkye1fPlyzZw5c8gLBwDY1zea\ncwEAIBzsLQYAsJwtwmWg/cns5sSJE3r88cc1f/58zZ8/X4WFhdq1a1e0ywqL3+/XxRdf7JjdEp55\n5hkVFhbqpptu0i9/+Uvt2bMn2iWF5PP5tGDBAt100026/vrrtWnTpmiX1Kf+9gr8z3/+o8LCQt14\n44264447bLPwua969+3bp3vvvVe33HKLCgsLtXz5cnONnR309//4tNLSUv3sZz8b5qr611+9X375\npX73u9+psLBQ+fn5uu2229TT0zPg9WwRLqH2J7Ojzz77TE8++aT+/Oc/64knnlBxcbGWLFmizz//\nPNqlDai0tFRnn312tMsIS01NjWpra/XUU0+prKxMt956q7744otolxXSypUrNWPGDJWVlWnr1q3a\ntm2bXn311WiXFaS/vQKPHz+uO++8U8uXL9dTTz2ladOmqaSkJEpV/q/+6n322WeVmJiobdu2aceO\nHRozZowt6pVC78coSe3t7aqoqOi1RjBaQtW7fPlyzZo1Szt27FBVVZXi4uKcES4D7U9mR/Hx8Sou\nLjb3WvvRj36kb33rW/rnP/8Z5cpCq6+vV3x8vPnYuN1t2rRJixcv1pgxp/6a5ufn6+qrr45yVaE1\nNzfrkksukXRqh4opU6boww8/jHJVwfrbK/C1115TTEyMfvjDH0qSbrjhBr344os6dOhQNMo09Vdv\nTk6ObrzxRvPzL37xC/3jH/8Y7vL61F/Np23atCmo9mjrr9733ntPfr9f11xzjdlWWlqqs84aeBVL\n1MNloP3J7CgpKUnXXXddUNvx48dt/4/2hg0bVFRUJCc8w3Hw4EF9+OGH8vv9uuWWW3TzzTfr6aef\njnZZA8rNzdXLL78s6dTf7Y8//lizZs2KclXB+tsrsL6+PujnMDU1VWPHjjV32YiW/uqdPXu20tLS\nzM/Hjh3ThAkThquskPqrWZL27t2rL774Qj/4wQ+GsaLQ+qt3165dSktL0wMPPKDCwkLdfvvtYf+y\nFPVwGWh/Mid48803NXnyZPM3Pjt67bXXdNFFF+m8886Ldilh+eSTTyRJzz//vB577DGVlpbq0Ucf\n1fPPPx/lykK7//771dLSoiuvvFIFBQX6/e9/b+4Mbnft7e29hkXGjx+v9nZnvKL31Vdf1U033RTt\nMgb0pz/9SXfddVe0ywjLvn37tHv3bk2dOlU7duzQddddpwULFph7RoYS9XA5rb/9yezu2LFj8nq9\neuCBB6JdSr8Mw9Cjjz6q22+/PdqlhO30i+UKCwsVExOj5ORkXXfddeaCXbsqKirS97//fb344ot6\n9tlntX79eu3evTvaZYWtrzkAJ/wsvvPOO2ptbdUtt9wS7VJCamhoUGxsrC666KJolxKW7u5uxcfH\n6/rrr5ckXXPNNYqNjdUrr7wy4PdGPVxSUlIkSYFAwGzr6uoy2+2upKREixYtUnZ2drRL6ddzzz2n\nH//4x0G9Q7s7vYP21/8enHvuufrss8+iVdKAfD6fdu/ebf4Dl5qaqjlz5jhiOE+SkpOT1dXVFdQW\nCARs/7O4b98+eb1elZaWKiYmJtrlhLRhwwaz1+KE0B4/fnyvocZzzz1XBw4cGPB7I9pbzArh7E9m\nV+vWrdP06dOVl5en7u5utbe369vf/na0y+plz549am5u1htvvCHDMNTe3q4//OEPmjRpkh5++OFo\nl9enKVOmaNy4cUHDo4cOHdK5554bxapCO378uCQFPY131llnOeb13t/73vdUU1Njfv7000917Ngx\n5eTkRLGq0A4dOqSVK1dq3bp1Sk5O1qeffqpzzjnHlk9EHj58WHv37tXKlSslnQrutrY2LViwQNdc\nc40KCwujXGFvfb1BONyfw6j3XJy6P9kjjzyiEydOKD8/X0eOHNHevXttO2SzZs0aPfXUU9q+fbue\neOIJpaSk6N5777VtsEhSbGys8vPzzWfuv/rqK73wwgtm99yOMjMzdd5552nnzp2STq0PePnll203\nod+fn/zkJzpx4oS5lqiyslJXXnmlbSbJz3TkyBEVFxeruLhYEyZM0OHDh7Vjxw7bPq4eHx+vmpoa\nbd++Xdu3b9e9996rc845R9u3b7dlsEjSnDlzdNZZZ6m2tlaS9Pbbb+vw4cOaPXv2gN9ri+1fnLY/\nmd/v11VXXWWOTxuGIZfLpTvvvFNFRUVRrq5/jY2N2rRpk15//XVdfPHFuvLKK7VgwYJol9Wvo0eP\n6r777tMHH3ygsWPHKi8vT4sWLYp2WSG99957WrdunQzD0JEjR3T55ZfrnnvuMR+ntoNQewU2Njbq\nvvvuU0xMjMaPH69169ZF/SV//dXr9Xr1yCOPmOed/jmsra3VpEmTolhx6P/HkrRt2za99NJLamho\n0E9+8hPdcccd+u53v2vLeuvr67VmzRqdffbZiomJ0cqVK8PqzdoiXAAAI4t9fp0CAIwYhAsAwHKE\nCwDAcoQLAMByhAsAwHKECwDAcoQLAMByhAsAwHKECwDAcv8PeZ7c8FmemfYAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f31616bb890>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.hist(M_ns.data, bins=50)\n",
    "plt.yscale('log')\n",
    "pass"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.23764238509276445"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "float(M_ns.nnz) / np.prod(M_ns.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Train the model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "scale = 0.03\n",
    "\n",
    "n_components = 100\n",
    "max_iter = 20\n",
    "n_jobs = 8\n",
    "lam_theta = lam_beta = 1e-5 * scale\n",
    "lam_gamma = 1e-5\n",
    "c0 = 1. * scale\n",
    "c1 = 10. * scale\n",
    "\n",
    "save_dir = os.path.join(DATA_DIR, 'ML20M_ns%d_scale%1.2E' % (k_ns, scale))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "reload(cofacto)\n",
    "coder = cofacto.CoFacto(n_components=n_components, max_iter=max_iter, batch_size=1000, init_std=0.01, n_jobs=n_jobs, \n",
    "                        random_state=98765, save_params=True, save_dir=save_dir, early_stopping=True, verbose=True, \n",
    "                        lam_theta=lam_theta, lam_beta=lam_beta, lam_gamma=lam_gamma, c0=c0, c1=c1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "ITERATION #0\n",
      "\tUpdating user factors: time=7.40\n",
      "\tUpdating item factors: time=6.99\n",
      "\tUpdating context factors: time=5.63\n",
      "\tUpdating bias terms: time=5.99\n",
      "\tValidation NDCG@k: 0.19163\n",
      "ITERATION #1\n",
      "\tUpdating user factors: time=5.65\n",
      "\tUpdating item factors: time=7.21\n",
      "\tUpdating context factors: time=5.19\n",
      "\tUpdating bias terms: time=5.94\n",
      "\tValidation NDCG@k: 0.26963\n",
      "ITERATION #2\n",
      "\tUpdating user factors: time=5.58\n",
      "\tUpdating item factors: time=7.58\n",
      "\tUpdating context factors: time=5.86\n",
      "\tUpdating bias terms: time=5.94\n",
      "\tValidation NDCG@k: 0.34485\n",
      "ITERATION #3\n",
      "\tUpdating user factors: time=5.93\n",
      "\tUpdating item factors: time=7.09\n",
      "\tUpdating context factors: time=5.57\n",
      "\tUpdating bias terms: time=5.90\n",
      "\tValidation NDCG@k: 0.35849\n",
      "ITERATION #4\n",
      "\tUpdating user factors: time=5.62\n",
      "\tUpdating item factors: time=7.35\n",
      "\tUpdating context factors: time=5.51\n",
      "\tUpdating bias terms: time=5.97\n",
      "\tValidation NDCG@k: 0.36298\n",
      "ITERATION #5\n",
      "\tUpdating user factors: time=5.52\n",
      "\tUpdating item factors: time=6.91\n",
      "\tUpdating context factors: time=5.61\n",
      "\tUpdating bias terms: time=5.96\n",
      "\tValidation NDCG@k: 0.36325\n",
      "ITERATION #6\n",
      "\tUpdating user factors: time=5.50\n",
      "\tUpdating item factors: time=7.12\n",
      "\tUpdating context factors: time=5.42\n",
      "\tUpdating bias terms: time=5.92\n",
      "\tValidation NDCG@k: 0.36150\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "CoFacto(batch_size=1000, dtype='float32', early_stopping=True, init_std=0.01,\n",
       "    max_iter=20, n_components=100, n_jobs=8, random_state=98765,\n",
       "    save_dir='/hdd2/dawen/data/ml-20m/pro/ML20M_ns1_scale3.00E-02',\n",
       "    save_params=True, verbose=True)"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "coder.fit(train_data, M_ns, vad_data=vad_data, batch_users=5000, k=100)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "test_data, _ = load_data(os.path.join(DATA_DIR, 'test.csv'))\n",
    "test_data.data = np.ones_like(test_data.data)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "n_params = len(glob.glob(os.path.join(save_dir, '*.npz')))\n",
    "\n",
    "params = np.load(os.path.join(save_dir, 'CoFacto_K%d_iter%d.npz' % (n_components, n_params - 1)))\n",
    "U, V = params['U'], params['V']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Test Recall@20: 0.1448\n",
      "Test Recall@50: 0.1765\n",
      "Test NDCG@100: 0.1724\n",
      "Test MAP@100: 0.0545\n"
     ]
    }
   ],
   "source": [
    "print 'Test Recall@20: %.4f' % rec_eval.recall_at_k(train_data, test_data, U, V, k=20, vad_data=vad_data)\n",
    "print 'Test Recall@50: %.4f' % rec_eval.recall_at_k(train_data, test_data, U, V, k=50, vad_data=vad_data)\n",
    "print 'Test NDCG@100: %.4f' % rec_eval.normalized_dcg_at_k(train_data, test_data, U, V, k=100, vad_data=vad_data)\n",
    "print 'Test MAP@100: %.4f' % rec_eval.map_at_k(train_data, test_data, U, V, k=100, vad_data=vad_data)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "np.savez('CoFactor_K100_ML20M.npz', U=U, V=V)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}
